C/Zaawansowane operacje matematyczne

Biblioteki

Biblioteka math.h
complex.h - liczby zespolone
liczby rzeczywiste
- float.h - liczny zmiennoprzecinkowe
- fenv.h - środowiski zmiennoprzecinkowe
liczby całkowite
- limits.h - charakterystyka typów liczb calkowitych
- inittypes.h - konwersja liczb całkowitych
- stdckdiniot.h - kontrola arytmetyki liczb całkowitych
atomic.h - typy atomowe
stdbit.h - marzędzia bitowe i bajtowe
stdbool.h - wartości i typy logiczne
tgmath.h - funkcje matematyczne niezależne od typu

Liczby

Liczby 0–F₁₆ w systemie dziesiętnym, ósemkowym i dwójkowym
0_hex	=	0_dec	=	0_oct	0	0	0	0
1_hex	=	1_dec	=	1_oct	0	0	0	1
2_hex	=	2_dec	=	2_oct	0	0	1	0
3_hex	=	3_dec	=	3_oct	0	0	1	1
4_hex	=	4_dec	=	4_oct	0	1	0	0
5_hex	=	5_dec	=	5_oct	0	1	0	1
6_hex	=	6_dec	=	6_oct	0	1	1	0
7_hex	=	7_dec	=	7_oct	0	1	1	1
8_hex	=	8_dec	=	10_oct	1	0	0	0
9_hex	=	9_dec	=	11_oct	1	0	0	1
A_hex	=	10_dec	=	12_oct	1	0	1	0
B_hex	=	11_dec	=	13_oct	1	0	1	1
C_hex	=	12_dec	=	14_oct	1	1	0	0
D_hex	=	13_dec	=	15_oct	1	1	0	1
E_hex	=	14_dec	=	16_oct	1	1	1	0
F_hex	=	15_dec	=	17_oct	1	1	1	1

Rodzaje liczb

całkowite
- int
- char
wymierne
- z pomocą typu mpq_t z biblioteki GMP ( GNU Multiple Precision Arithmetic Library )
- z użyciem struktur ^[1]^[2]
zmiennoprzecinkowe
- float
- double
- wartości specjalne
zespolone
kwaterniony

Dodatkowo rodzaj liczby definiują :

specyfikatory
podstawa systemu liczbowego
- binarne
- ósemkowe (oktalne)
- dziesiętne
- szesnastkowe

Jak widać nie ma :

Liczb niewymiernych ( chyba że korzystamy z obliczeń symbolicznych, np. program Maxima CAS )

liczby binarne

Sposoby :

itoa
snprintf ^[3]
Binarna stała całkowita z użyciem rozszerzenia gcc
za pomoc operacji na bitach, np. shift

Binarna stała całkowita z użyciem rozszerzenia gcc ^[4]

prefix ‘0b’ lub ‘0B’
sekwencja cyfr ‘0’ lub ‘1’
suffix ‘L’ lub ‘UL’

Następujące zapisy są identyczne:

     i =       42; // decimal
     i =     0x2a; // hexadecimal
     i =      052; // octal
     i = 0b101010; // binary

int i = 1 << 9; /* binary 1 followed by 9 binary 0's */ ^[5]

Liczby szestnastkowe ( hexadecymalne)

Liczba szesnastkowa:

jest poprzedzana przez przedrostek „0x” lub „0X”
jest wartością całkowitą i można ją przechowywać w typie całkowitym typów danych (char, short lub int)

Przykłady:

„0x64” odpowiada 100 w systemie dziesiętnym
wartość szesnastkowa 3E8 jest reprezentowana w c jako 0x3E8

Każdy bajt ( 8 bitów) to 2-cyfrowa liczba szesnastkowa

Maksymalna liczba całkowita bez znaku zajmująca n bajtów:

1 bajt = 8 bitów: $FF_{16}=11111111_{2}=255_{10}=2^{8}-1$
2 bajty = 16 bitów $FFFF_{16}=65535_{10}=2^{16}-1$
4 bajty = 32 bity $FFFFFFFF_{16}=4294967295_{10}=2^{32}-1$
8 bajtów = 64 bity $FFFFFFFFFFFFFFFF_{16}=18446744073709551615_{10}=2^{64}-1$

Prezentacja liczb rzeczywistych w pamięci komputera

Liczby rzeczywisty

W wielu książkach nie ma w ogóle tego tematu. Być może ten temat może wydać Ci się niepotrzebnym, lecz dzięki niemu zrozumiesz, jak komputer radzi sobie z przecinkiem ^[6]oraz dlaczego niektóre obliczenia dają niezbyt dokładne wyniki.^[7]
Na początek trochę teorii. Do przechowywania liczb rzeczywistych przeznaczone są 3 typy: float, double oraz long double. Zajmują one odpowiednio 32, 64 oraz 80 bitów. Wiemy też, że komputer nie ma fizycznej możliwości zapisania przecinka. Spróbujmy teraz zapisać jakąś liczbę wymierną w formie liczb binarnych. Nasza liczba to powiedzmy 4.25. Spróbujmy ją rozbić na sumę potęg dwójki: 4 = 1*2² + 0*2¹+0*2⁰. Dobra - rozpisaliśmy liczbę 4, ale co z częścią dziesiętną? Skorzystajmy z zasad matematyki - 0.25 = 2^-2. Zatem nasza liczba powinna wyglądać tak:

100.01

Ponieważ komputer nie jest w stanie przechować pozycji przecinka, ktoś wpadł na prosty ale sprytny pomysł ustawienia przecinka jak najbliżej początku liczby i tylko mnożenia jej przez odpowiednią potęgę dwójki. Taki sposób przechowywania liczb nazywamy zmiennoprzecinkowym, a proces przekształcania naszej liczby z postaci czytelnej przez człowieka na format zmiennoprzecinkowy nazywamy normalizacją. Wróćmy do naszej liczby - 4.25. W postaci binarnej wygląda ona tak: 100.01, natomiast po normalizacji będzie wyglądała tak: 1.0001*2². W ten sposób w pamięci komputera znajdą się dwie informacje: liczba zakodowana w pamięci z "wirtualnym" przecinkiem oraz numer potęgi dwójki. Te dwie informacje wystarczają do przechowania wartości liczby. Jednak pojawia się inny problem - co się stanie, jeśli np. będziemy chcieli przełożyć liczbę typu ${\frac {1}{3}}$ ? Otóż tutaj wychodzą na wierzch pewne niedociągnięcia komputera w dziedzinie samej matematyki. 1/3 daje w rozwinięciu dziesiętnym 0.(3). Jak zatem zapisać taką liczbę? Otóż nie możemy przechować całego jej rozwinięcia (wynika to z ograniczeń typu danych - ma on niestety skończoną liczbę bitów). Dlatego przechowuje się tylko pewne przybliżenie liczby. Jest ono tym bardziej dokładne im dany typ ma więcej bitów. Zatem do obliczeń wymagających dokładnych danych powinniśmy użyć typu double lub long double. Na szczęście w większości przeciętnych programów tego typu problemy zwykle nie występują. A ponieważ początkujący programista nie odpowiada za tworzenie programów sterujących np. lotem statku kosmicznego, więc drobne przekłamania na odległych miejscach po przecinku nie stanowią większego problemu.

Program wyświetlający wewnętrzną reprezentację liczby podwójnej precyzji:^[8]^[9]

#include <stdio.h>

int main(void)
{
  double a = 1.0 / 3;
  size_t i;
  size_t iMax= sizeof a;

  printf("bytes are numbered from 0 to %x\n", (unsigned)iMax -1);

  for (i = 0; i < iMax ; ++i) 
  {
    printf("byte number %u  =  %x\n", (unsigned)i, ((unsigned char *)&a)[i]);
  }


  printf("hex memory representation  of %f is : \n", a);
  for (i = iMax; i>0 ; --i) 
  {
    printf("%x", ((unsigned char *)&a)[i-1]);
  }
  printf(" \n");
  return 0;
}

Daje wynik:

bytes are numbered from 0 to 7
byte number 0  =  55
byte number 1  =  55
byte number 2  =  55
byte number 3  =  55
byte number 4  =  55
byte number 5  =  55
byte number 6  =  d5
byte number 7  =  3f

hex memory representation  of 0.333333 is : 
3fd5555555555555

0 01111111101 0101010101010101010101010101010101010101010101010101₂ = 3FD5 5555 5555 5555₁₆ ≙ +2⁻² × (1 + 2⁻² + 2⁻⁴ + ... + 2⁻⁵²) ≈ ¹/₃

Given the hexadecimal representation 3FD5 5555 5555 5555₁₆,
  Sign = 0
  Exponent = 3FD₁₆ = 1021
  Exponent Bias = 1023 (constant value; see above)
  Fraction = 5 5555 5555 5555₁₆
  Value = 2^{(Exponent − Exponent Bias)} × 1.Fraction – Note that Fraction must not be converted to decimal here
        = 2⁻² × (15 5555 5555 5555₁₆ × 2⁻⁵²)
        = 2⁻⁵⁴ × 15 5555 5555 5555₁₆
        = 0.333333333333333314829616256247390992939472198486328125
        ≈ 1/3

Zobacz również:

dumpfp^[10]

Liczby całkowite

Liczba całkowitej bez znaku (ang. unsigned char lub uint8_t ) zajmującej 1 bajt, czyli 8 bitów.

Ta tabela ilustruje wartość

 $149_{10}=10010101_{2}=95_{16}=225_{8}$

Oznaczenia

MSB = najbardziej znaczący bit
LSB = najmniej znaczący bit

Liczba dziesiętna bez znaku 149 w postaci binarnej z zaznaczonym LSB

Liczba dziesiętna bez znaku 149 w postaci binarnej z zaznaczonym MSB

$149=\sum _{i=0}^{7}b_{i}\cdot 2^{i}=1*2^{0}+1*2^{2}+1*2^{4}+1*2^{7}=1+4+16+128=149$

Bit position i	7	6	5	4	3	2	1	0
Binary digit $b_{i}$	1	0	0	1	0	1	0	1
Bit weight = ( 2ⁱ )	2⁷	2⁶	2⁵	2⁴	2³	2²	2¹	2⁰
Bit position label	MSB							LSB

Liczba całkowita bez znaku zajmująca 2 bajty = 16 bitów ( uint16_t )

Wartość

dziesiętna (decymalna) = 769
dwójkowa (binarna) = 1100000001
hexadecymalna = 0x0301

$769_{10}=0301_{16}=1100000001_{2}$

 +---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+
 | 0 | 0 | 0 | 0 | 0 | 0 | 1 | 1 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 1 |  
 +---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+

Zadania:

konwersja dwóch wartości 8-bitowych na jedną wartość 16-bitową w C^[11]
dzielenie całkowitej liczby 16 bitowej bez znaku na 2 bajty^[12]

Konwersje liczb

/*

gcc b.c -Wall -Wextra


*/


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


/*
https://stackoverflow.com/questions/5488377/converting-an-integer-to-binary-in-c


If you want to transform a number into another number (not number to string of characters), and you can do with a small range (0 to 1023 for implementations with 32-bit integers), you don't need to add char* to the solution. pmg

*/
unsigned unsigned_int_to_bin(unsigned int k) {
    if (k == 0) return 0;
    if (k == 1) return 1;                       /* optional */
    return (k % 2) + 10 * unsigned_int_to_bin(k / 2);
}



int main(void){


	unsigned int u = 1019; 
	if (u>1023) {fprintf(stderr, "input number decimal u = %d is too big; close the program ! \n", u); return 1;}
	
	unsigned int b = unsigned_int_to_bin(u);
	fprintf(stdout, "decimal u = %d as a binary unsigned number = %d\n", u, b);
	



	return 0;

}

Wynik:

gcc b.c -Wall -Wextra
./a.out
decimal u = 1019 as a binary unsigned number = 1111111011

Konwersja zapisu matematycznego na program komputerowy w języku C

sumowanie

W zapisie matematycznym do przedstawiania w zwarty sposób sumowania wielu podobnych wyrazów ( serii) stosuje się symbol sumowania $\Sigma$ , wywodzący się z wielkiej greckiej litery sigma. Jest on zdefiniowany jako:

\sum _{i=m}^{n}x_{i}=x_{m}+x_{m+1}+x_{m+2}+\ldots +x_{n-1}+x_{n},

gdzie

$i$ oznacza indeks sumowania
$x_{i}$ to zmienna przedstawiająca każdy kolejny wyraz w szeregu
$m$ to dolna granica sumowania
$n$ to górna granica sumowania.

Wyrażenie „i = m” pod symbolem sumowania oznacza, że indeks $i$ rozpoczyna się od wartości $m.$ Następnie dla każdego kolejnego wyrazu indeks $i$ jest zwiększany o 1, aż $i$ osiągnie wartość $n$ (tj. $i=n$ ), który jest ostatnim wyrazem sumowania.

#include <stdio.h>

int summation(const int m,  const int n )
{
   int s = 0;
   for(int i=m; i<=n; ++i)
   {
     s += i;
   }
   return s;
}



int main()
{
   int sum;
   int m = 1; // lower index
   int n = 100; // upper index
   
   sum = summation(m,n);
   
   printf("sum of integer numbers from %d to %d is %d\n", m, n, sum);
   return 0; 
 }

 gcc s.c -Wall -Wextra
 ./a.out
 sum of integer numbers from 1 to 100 is 5050

Suma sum:

$\sum _{i=1}^{4}\sum _{j=1}^{2}ij^{2}=50$

Sumy są obliczane od prawej do lewej. Oznacza to, że sigma po lewej stronie jest najwyżej zagnieżdżoną pętlą: ^[13]

    
sum = 0;
for (i =1; i<=4; i++)
    for (j=1; j <=2; j++)
         sum += i * j* j;

produkt ( iloczyn)

W zapisie matematycznym do przedstawiania w zwarty sposób iloczynu wielu podobnych wyrazów ( serii) stosuje się symbol iloczynu $\Pi$ , wywodzący się z wielkiej greckiej litery pi. Jest on zdefiniowany jako:

\prod _{i=1}^{n}x_{i}=x_{1}\cdot x_{2}\cdot x_{3}\cdot \ldots \cdot x_{n}

#include <stdio.h>

int product(const int m,  const int n )
{
   int p = 1;
   for(int i=m; i<=n; ++i)
   {
     p *= i;
   }
   return p;
}



int main()
{
   int p;
   int m = 1; // lower index
   int n = 10; // upper index
   
    p = product(m,n);
   
   printf("product of integer numbers from %d to %d is %d\n", m, n, p);
   
}

Dla zakresu [1,10] program działa poprawnie

gcc p.c -Wall -Wextra
./a.out
product of integer numbers from 1 to 10 is 3628800

Dla zakresu [1,100] program daje dziwny wynik:

gcc p.c -Wall -Wextra
./a.out
product of integer numbers from 1 to 100 is 0

Wartość iloczynu możemy obliczyć za pomoca wolfram alfa.

Product[k, {k, 1, 100}] = 93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000

Liczba ma 159 cyfr, to więcej niż zakres dostępnych standardowych liczb całkowitych. Przekroczenie zakresu powoduje błędny wynik.

Jeśli dolny zakres wynosi 1 to wynik jest równy silni ( ang factorial). Wartość silnii rośnie szybciej niż wzrost funkcji wykładniczej.

n!=\prod _{i=1}^{n}i.

Sprawdzenie czy nie będzie przekroczenia zakresu liczb powinniśmy wykonać przed wykonaniem działania:

#include <stdio.h>
#include <limits.h> // INT_MAX
/*
http://c-faq.com/misc/sd26.html
functions for ``careful'' multiplication.

*/

int chkmul(int a, int b)
{
	int sign = 1;
	if(a == 0 || b == 0) return 0;
	if(a < 0) { a = -a; sign = -sign; }
	if(b < 0) { b = -b; sign = -sign; }
	if(INT_MAX / b < a) {
		fprintf(stdout,"int overflow\t");
		return (sign > 0) ? INT_MAX : INT_MIN;
	}
	
	
	return sign*a*b;
}


int product(const int m,  const int n )
{
   int p = 1;
   for(int i=m; i<=n; ++i)
   {
     fprintf(stdout,"i = %d \t", i);
     p = chkmul(p,i);
     fprintf(stdout,"p = %d \n", p);
     
   }
   return p;
}



int main()
{
   int p;
   int m = 1; // lower index
   int n = 100; // upper index
   
    p = product(m,n);
   
   printf("product of integer numbers from %d to %d is %d\n", m, n, p);
   
}

gcc p.c -Wall -Wextra
./a.out
i = 1 	p = 1 
i = 2 	p = 2 
i = 3 	p = 6 
i = 4 	p = 24 
i = 5 	p = 120 
i = 6 	p = 720 
i = 7 	p = 5040 
i = 8 	p = 40320 
i = 9 	p = 362880 
i = 10 	p = 3628800 
i = 11 	p = 39916800 
i = 12 	p = 479001600 
i = 13 	int overflow	p = 2147483647 
i = 14 	int overflow	p = 2147483647 
i = 15 	int overflow	p = 2147483647 
i = 16 	int overflow	p = 2147483647 
i = 17 	int overflow	p = 2147483647 
i = 18 	int overflow	p = 2147483647 
i = 19 	int overflow	p = 2147483647 
i = 20 	int overflow	p = 2147483647 
i = 21 	int overflow	p = 2147483647 
i = 22 	int overflow	p = 2147483647 
i = 23 	int overflow	p = 2147483647 
i = 24 	int overflow	p = 2147483647 
i = 25 	int overflow	p = 2147483647 
i = 26 	int overflow	p = 2147483647 
i = 27 	int overflow	p = 2147483647 
i = 28 	int overflow	p = 2147483647 
i = 29 	int overflow	p = 2147483647 
i = 30 	int overflow	p = 2147483647 
i = 31 	int overflow	p = 2147483647 
i = 32 	int overflow	p = 2147483647 
i = 33 	int overflow	p = 2147483647 
i = 34 	int overflow	p = 2147483647 
i = 35 	int overflow	p = 2147483647 
i = 36 	int overflow	p = 2147483647 
i = 37 	int overflow	p = 2147483647 
i = 38 	int overflow	p = 2147483647 
i = 39 	int overflow	p = 2147483647 
i = 40 	int overflow	p = 2147483647 
i = 41 	int overflow	p = 2147483647 
i = 42 	int overflow	p = 2147483647 
i = 43 	int overflow	p = 2147483647 
i = 44 	int overflow	p = 2147483647 
i = 45 	int overflow	p = 2147483647 
i = 46 	int overflow	p = 2147483647 
i = 47 	int overflow	p = 2147483647 
i = 48 	int overflow	p = 2147483647 
i = 49 	int overflow	p = 2147483647 
i = 50 	int overflow	p = 2147483647 
i = 51 	int overflow	p = 2147483647 
i = 52 	int overflow	p = 2147483647 
i = 53 	int overflow	p = 2147483647 
i = 54 	int overflow	p = 2147483647 
i = 55 	int overflow	p = 2147483647 
i = 56 	int overflow	p = 2147483647 
i = 57 	int overflow	p = 2147483647 
i = 58 	int overflow	p = 2147483647 
i = 59 	int overflow	p = 2147483647 
i = 60 	int overflow	p = 2147483647 
i = 61 	int overflow	p = 2147483647 
i = 62 	int overflow	p = 2147483647 
i = 63 	int overflow	p = 2147483647 
i = 64 	int overflow	p = 2147483647 
i = 65 	int overflow	p = 2147483647 
i = 66 	int overflow	p = 2147483647 
i = 67 	int overflow	p = 2147483647 
i = 68 	int overflow	p = 2147483647 
i = 69 	int overflow	p = 2147483647 
i = 70 	int overflow	p = 2147483647 
i = 71 	int overflow	p = 2147483647 
i = 72 	int overflow	p = 2147483647 
i = 73 	int overflow	p = 2147483647 
i = 74 	int overflow	p = 2147483647 
i = 75 	int overflow	p = 2147483647 
i = 76 	int overflow	p = 2147483647 
i = 77 	int overflow	p = 2147483647 
i = 78 	int overflow	p = 2147483647 
i = 79 	int overflow	p = 2147483647 
i = 80 	int overflow	p = 2147483647 
i = 81 	int overflow	p = 2147483647 
i = 82 	int overflow	p = 2147483647 
i = 83 	int overflow	p = 2147483647 
i = 84 	int overflow	p = 2147483647 
i = 85 	int overflow	p = 2147483647 
i = 86 	int overflow	p = 2147483647 
i = 87 	int overflow	p = 2147483647 
i = 88 	int overflow	p = 2147483647 
i = 89 	int overflow	p = 2147483647 
i = 90 	int overflow	p = 2147483647 
i = 91 	int overflow	p = 2147483647 
i = 92 	int overflow	p = 2147483647 
i = 93 	int overflow	p = 2147483647 
i = 94 	int overflow	p = 2147483647 
i = 95 	int overflow	p = 2147483647 
i = 96 	int overflow	p = 2147483647 
i = 97 	int overflow	p = 2147483647 
i = 98 	int overflow	p = 2147483647 
i = 99 	int overflow	p = 2147483647 
i = 100 	int overflow	p = 2147483647 
product of integer numbers from 1 to 100 is 2147483647

Widać że każdy produkt powyżej górnego zakresu=12 jest błędny. Rozwiązaniem jest:

użycie bibliotek o dowolnej precyzji ( GMP )
użycie typu zmiennoprzecinkowego (double )

#include <stdio.h>

double double_product(const int m,  const int n )
{
   double p = 1.0;
   for(int i=m; i<=n; ++i)
   {
     p*=i;
     
   }
   return p;
}


int main()
{
   double p;
   int m = 1; // lower index
   int n = 100; // upper index
   
    p = double_product(m,n);
   
   printf("product of integer numbers from %d to %d is %.16e\n", m, n, p);
   
}

gcc d.c -Wall -Wextra
./a.out
product of integer numbers from 1 to 100 is 9.3326215443944102e+157

Jak widać pierwsze 17 cyfr znaczących się zgadza. Błąd wynosi około 10^140 czyli jest mniejszy niż 1%

zaawansowane algorytmy

Jak tworzyć matematyczne oprogramowanie w C?^[14]
drukowanie liczb algebraicznych^[15]

Obliczenia numeryczne

Obliczenia numeryczne^[16] ^[17]są to obliczenia na liczbach. Ich przeciwieństwem są obliczenia symboliczne wykonywane na symbolach ( zobacz Maxima CAS ) ^[18]

Epsilon maszynowy

Epsilon maszynowy jest wartością określającą precyzję obliczeń numerycznych wykonywanych na liczbach zmiennoprzecinkowych.^[21]

Jest to najmniejsza liczba nieujemna, której dodanie do jedności daje wynik nierówny 1. Innymi słowy, jest to najmniejszy ε, dla którego następujący warunek jest uznawany za niespełniony (przyjmuje wartość fałsz): 1 + ε = 1

Im mniejsza wartość epsilona maszynowego, tym większa jest względna precyzja obliczeń.

Obliczmy epsilon dla liczb podwójnej precyzji :

/* 
http://en.wikipedia.org/wiki/Machine_epsilon
The following C program does not actually determine the machine epsilon;
rather, it determines a number within a factor of two (one order of magnitude) 
of the true machine epsilon, using a linear search.

---

The difference between 1 and the least value greater than 1 that is representable in the given floating-point type, b1-p.
-------------------------------
http://stackoverflow.com/questions/1566198/how-to-portably-get-dbl-epsilon-in-c-c

gcc m.c -Wall

./a.out
*/


#include <stdio.h>
#include <float.h> // DBL_EPSILON
 
 int main()
 {
    double epsilon = 1.0;
 
    printf( "epsilon;  1 + epsilon\n" );
    
    do 
     {
       printf( "%G\t%.20f\n", epsilon, (1.0 + epsilon) );
       epsilon /= 2.0f;
     }
    // If next epsilon yields 1, then break
    while ((1.0 + (epsilon/2.0)) != 1.0); // 

    // because current epsilon is the calculated machine epsilon.
    printf( "\nCalculated Machine epsilon: %G\n", epsilon );
    
    
    //check value from float.h , Steve Jessop
    if ((1.0 + DBL_EPSILON) != 1.0 
        && 
       (1.0 + DBL_EPSILON/2) == 1.0)
     printf("DBL_EPSILON = %g \n", DBL_EPSILON);
     else printf("DBL_EPSILON is not good !!! \n");


    return 0;
 }

Wynik programu :

epsilon;  1 + epsilon
1	2.00000000000000000000
0.5	1.50000000000000000000
0.25	1.25000000000000000000
0.125	1.12500000000000000000
0.0625	1.06250000000000000000
0.03125	1.03125000000000000000
0.015625	1.01562500000000000000
0.0078125	1.00781250000000000000
0.00390625	1.00390625000000000000
0.00195312	1.00195312500000000000
0.000976562	1.00097656250000000000
0.000488281	1.00048828125000000000
0.000244141	1.00024414062500000000
0.00012207	1.00012207031250000000
6.10352E-05	1.00006103515625000000
3.05176E-05	1.00003051757812500000
1.52588E-05	1.00001525878906250000
7.62939E-06	1.00000762939453125000
3.8147E-06	1.00000381469726562500
1.90735E-06	1.00000190734863281250
9.53674E-07	1.00000095367431640625
4.76837E-07	1.00000047683715820312
2.38419E-07	1.00000023841857910156
1.19209E-07	1.00000011920928955078
5.96046E-08	1.00000005960464477539
2.98023E-08	1.00000002980232238770
1.49012E-08	1.00000001490116119385
7.45058E-09	1.00000000745058059692
3.72529E-09	1.00000000372529029846
1.86265E-09	1.00000000186264514923
9.31323E-10	1.00000000093132257462
4.65661E-10	1.00000000046566128731
2.32831E-10	1.00000000023283064365
1.16415E-10	1.00000000011641532183
5.82077E-11	1.00000000005820766091
2.91038E-11	1.00000000002910383046
1.45519E-11	1.00000000001455191523
7.27596E-12	1.00000000000727595761
3.63798E-12	1.00000000000363797881
1.81899E-12	1.00000000000181898940
9.09495E-13	1.00000000000090949470
4.54747E-13	1.00000000000045474735
2.27374E-13	1.00000000000022737368
1.13687E-13	1.00000000000011368684
5.68434E-14	1.00000000000005684342
2.84217E-14	1.00000000000002842171
1.42109E-14	1.00000000000001421085
7.10543E-15	1.00000000000000710543
3.55271E-15	1.00000000000000355271
1.77636E-15	1.00000000000000177636
8.88178E-16	1.00000000000000088818
4.44089E-16	1.00000000000000044409

Calculated Machine epsilon: 2.22045E-16
DBL_EPSILON = 2.22045e-16

Obliczmy epsilon dla liczb pojedynczej precyzji :

#include <stdio.h>
 
 int main()
 {
    float epsilon = 1.0f;
 
    printf( "epsilon;  1 + epsilon\n" );
    
    do 
     {
       printf( "%G\t%.20f\n", epsilon, (1.0f + epsilon) );
       epsilon /= 2.0f;
     }
    // If next epsilon yields 1, then break
    while ((float)(1.0 + (epsilon/2.0)) != 1.0); // 

    // because current epsilon is the machine epsilon.
    printf( "\nCalculated Machine epsilon: %G\n", epsilon );
    return 0;
 }

Wynik :

Calculated Machine epsilon: 1.19209E-07

Obliczmy epsilon dla liczb long double :

#include <stdio.h>
 
 int main()
 {
    long double epsilon = 1.0;
 
    printf( "epsilon;  1 + epsilon\n" );
    
    do 
     {
       printf( "%LG \t %.25LG \n", epsilon, (1.0 + epsilon) );
       epsilon /= 2.0;
     }
    // If next epsilon yields 1, then break
    while ((1.0 + (epsilon/2.0)) != 1.0); // 

    // because current epsilon is the machine epsilon.
    printf( "\n Calculated Machine epsilon: %LG\n", epsilon );
    return 0;
 }

Wynik :

 Calculated Machine epsilon: 1.0842E-19

Porównaj

wersja angielska

wyjątki

Wyjatki

programowe
sprzętowe

Przykłady

A floating-point exception (FPE) = Wyjątek zmiennoprzecinkowy nie jest wyjątkiem programowym. Powstaje na poziomie mikroprocesora lub ISA. FPE może spowodować wysłąnie sygnału o nazwie SIGFPE, z którym można sobie poradzić, ale nie za pomocą C^[22]

Wyjątki zmiennoprzecinkowe są kontrolowane przez kod biblioteki w C99, a nie przez flagi kompilatora^[23]

    
#include <fenv.h>
#include <math.h>
#include <stdio.h>

#define PRINTEXC(ex, val) printf(#ex ": %s\n", (val & ex) ? "set" : "unset");

double foo(double a, double b) { return sin(a) / b; }

int main()
{
    int e;
    double x;

    feclearexcept(FE_ALL_EXCEPT);

    x = foo(1.2, 3.1);

    e = fetestexcept(FE_ALL_EXCEPT);
    PRINTEXC(FE_DIVBYZERO, e);
    PRINTEXC(FE_INEXACT, e);
    PRINTEXC(FE_INVALID, e);
    PRINTEXC(FE_OVERFLOW, e);
    PRINTEXC(FE_UNDERFLOW, e);

    putchar('\n');

    feclearexcept(FE_ALL_EXCEPT);

    x += foo(1.2, 0.0);

    e = fetestexcept(FE_ALL_EXCEPT);
    PRINTEXC(FE_DIVBYZERO, e);
    PRINTEXC(FE_INEXACT, e);
    PRINTEXC(FE_INVALID, e);
    PRINTEXC(FE_OVERFLOW, e);
    PRINTEXC(FE_UNDERFLOW, e);
    return lrint(x);
}

Wynik:

    

FE_DIVBYZERO: unset
FE_INEXACT: set
FE_INVALID: unset
FE_OVERFLOW: unset
FE_UNDERFLOW: unset

FE_DIVBYZERO: set
FE_INEXACT: set
FE_INVALID: unset
FE_OVERFLOW: unset
FE_UNDERFLOW: unset

Jest też możliwe użycie sygnałów ( zobacz fenv.h):

    
#pragma STDC FENV_ACCESS on

#define _GNU_SOURCE
#include <fenv.h>

int main()
{
#ifdef FE_NOMASK_ENV
    fesetenv(FE_NOMASK_ENV);
#endif

    // ...
}

Limity dla obliczeń

zmiennoprzecinkowych

Definicje W pliku float.h są zdefiniowane stałe :^[24]

DBL_MIN , czyli najmniejszą dodatnia liczbą typu double uznawaną przez maszynę za różną od zera ^[25]
DBL_MAX, czyli największa dodatnia liczbą typu double, która może być używana przez komputer

W pliku math.h są zdefiniowane :

NAN

// gcc -lm -Wall l.c
#include <stdio.h>
#include <math.h> // infinity, nan
#include <float.h>//DBL_MIN

int main(void)
{
  
  printf("DBL_MIN = %g \n", DBL_MIN);
  printf("DBL_MAX = %g \n", DBL_MAX);
  printf("INFINITY = %g \n",  INFINITY);
#ifdef NAN 
  printf("NAN= %g \n", NAN );
#endif
  return 0;
}

Wynik działania :

DBL_MIN = 2.22507e-308 
DBL_MAX = 1.79769e+308 
INFINITY = inf 
NAN= nan

całkowitych

/*

gcc l.c -lm -Wall
./a.out


http://stackoverflow.com/questions/29592898/do-long-long-and-long-have-same-range-in-c-in-64-bit-machine
*/
#include <stdio.h>
#include <math.h> // M_PI; needs -lm also 
#include <limits.h> // INT_MAX, http://pubs.opengroup.org/onlinepubs/009695399/basedefs/limits.h.html





int main(){

double lMax;


 lMax = log2(INT_MAX);
 printf("INT_MAX \t= %25d ; lMax = log2(INT_MAX) \t= %.0f \n",INT_MAX,  lMax); 

 lMax = log2(UINT_MAX);
 printf("UINT_MAX \t= %25u ; lMax = log2(UINT_MAX) \t= %.0f \n", UINT_MAX,  lMax); 

 lMax = log2(LONG_MAX);
 printf("LONG_MAX \t= %25ld ; lMax = log2(LONG_MAX) \t= %.0f \n",LONG_MAX,  lMax); 


 lMax = log2(ULONG_MAX);
 printf("ULONG_MAX \t= %25lu ; lMax = log2(ULONG_MAX) \t= %.0f \n",ULONG_MAX,  lMax); 

 lMax = log2(LLONG_MAX);
 printf("LLONG_MAX \t= %25lld ; lMax = log2(LLONG_MAX) \t= %.0f \n",LLONG_MAX, lMax); 

 lMax = log2(ULLONG_MAX);
 printf("ULLONG_MAX \t= %25llu ; lMax = log2(ULLONG_MAX) \t= %.0f \n",ULLONG_MAX, lMax); 





return 0;
}

Wynik :

INT_MAX 	=                2147483647 ; lMax = log2(INT_MAX) 	= 31 
UINT_MAX 	=                4294967295 ; lMax = log2(UINT_MAX) 	= 32 
LONG_MAX 	=       9223372036854775807 ; lMax = log2(LONG_MAX) 	= 63 
ULONG_MAX 	=      18446744073709551615 ; lMax = log2(ULONG_MAX) 	= 64 
LLONG_MAX 	=       9223372036854775807 ; lMax = log2(LLONG_MAX) 	= 63 
ULLONG_MAX 	=      18446744073709551615 ; lMax = log2(ULLONG_MAX) 	= 64

Dla typów o stałej szerokości ( ang. Fixed width integer types ):

/*

gcc h.c -Wall -Wextra
./a.out

printf
------------------------------------ 
X ( specifier character): 
	input  = argument of type ‘unsigned int’
	output =  Unsigned hexadecimal integer (uppercase) without 0X prefix



# (  flags sub-specifier)	Used with X specifiers  : the value is preceeded with 0X respectively for values different than zero


*/
#include <stdio.h>
#include <stdint.h>
#include <inttypes.h> // PRIu32 = all those nifty new format specifiers for the intN_t types and their brethren


int main(void)
{
	
	uint8_t  a = 0XFF; 		// 1 byte  =  8 bits
	uint16_t b = 0XFFFF;  		// 2 bytes = 16 bits
	uint32_t c = 0XFFFFFFFF;	// 4 bytes = 32 bits
	uint64_t d = 0XFFFFFFFFFFFFFFFF;// 8 bytes = 64 bits
	
	int bytes = 1;
	int bits = 8* bytes;
	
	
	printf("uint%d_t \ta ( %d byte  =  %d bits)\t = %#X (hexadecimal)\t\t\t= %d \t\t\t = 2^%d - 1  (decimal) \n", bits, bytes, bits, a, a, bits);
	bytes = 2; bits = 8 * bytes;
	printf("uint%d_t\tb ( %d byte  =  %d bits)\t = %#X (hexadecimal)\t\t\t= %d\t\t\t = 2^%d - 1  (decimal) \n", bits, bytes, bits, b, b, bits);
	bytes = 4; bits = 8 * bytes;
	printf("uint%d_t\tc ( %d byte  =  %d bits)\t = %#X (hexadecimal)\t\t= %"PRIu32" \t\t = 2^%d - 1  (decimal) \n", bits, bytes, bits, c, c, bits);
	bytes = 8; bits = 8 * bytes;
	printf("uint%d_t\tc ( %d byte  =  %d bits)\t = %#lX (hexadecimal)\t= %"PRIu64" \t = 2^%d - 1  (decimal) \n", bits, bytes, bits, d, d, bits);
	
	
	return 0;
}

Wynik:

uint8_t 	a ( 1 byte  =  8 bits)	 = 0XFF (hexadecimal)			= 255 			 = 2^8 - 1  (decimal) 
uint16_t	b ( 2 byte  =  16 bits)	 = 0XFFFF (hexadecimal)			= 65535			 = 2^16 - 1  (decimal) 
uint32_t	c ( 4 byte  =  32 bits)	 = 0XFFFFFFFF (hexadecimal)		= 4294967295 		 = 2^32 - 1  (decimal) 
uint64_t	c ( 8 byte  =  64 bits)	 = 0XFFFFFFFFFFFFFFFF (hexadecimal)	= 18446744073709551615 	 = 2^64 - 1  (decimal)

Limity dla typów całkowitych o stałej szerokości
Typ	Znak	Bits	Bytes	Minimum	Maximum
`int8_t`	Signed	8	1	−2⁷ = −128	2⁷ − 1 = 127
`uint8_t`	Unsigned	8	1	0	2⁸ − 1 = 255
`int16_t`	Signed	16	2	−2¹⁵ = −32,768	2¹⁵ − 1 = 32,767
`uint16_t`	Unsigned	16	2	0	2¹⁶ − 1 = 65,535
`int32_t`	Signed	32	4	−2³¹ = −2,147,483,648	2³¹ − 1 = 2,147,483,647
`uint32_t`	Unsigned	32	4	0	2³² − 1 = 4,294,967,295
`int64_t`	Signed	64	8	−2⁶³ = −9,223,372,036,854,775,808	2⁶³ − 1 = 9,223,372,036,854,775,807
`uint64_t`	Unsigned	64	8	0	2⁶⁴ − 1 = 18,446,744,073,709,551,615

Przekroczenie zakresu liczb całkowitych

Przekroczenie zakresu liczb całkowitych ( ang. integer overflow ) ^[26] może dotyczyć liczb całkowitych :^[27]

bez znaku ( " Unsigned integers are defined to wrap around. " )
ze znakiem ( powoduje zachowanie niezdefiniowane - może to powodować Złe Rzeczy czyli zagrożenie bezpieczeństwa komputera ^[28] )
- abs(INT_MIN)^[29]

#include <stdio.h>

/* 

a signed integer overflow is undefined behaviour in C 
check b^i 

to compile : 
  gcc i.c -Wall
to run : 

  ./a.out

*/

int main() {

int i;
int b=2; // base 
int p=1; // p = b^i

for ( i=0 ; i<40; i++){
     
     printf(" b^i = %d ^ %d = %d \n", b, i, p);
     p *= b;
     
}


return 0; 
}

Program kompiluje się i uruchamia bez komunikatów o błędach ale wynik nie jest taki jak naiwnie moglibyśmy się spodziewać :


 b^i = 2 ^ 0 = 1 
 b^i = 2 ^ 1 = 2 
 b^i = 2 ^ 2 = 4 
 b^i = 2 ^ 3 = 8 
 b^i = 2 ^ 4 = 16 
 b^i = 2 ^ 5 = 32 
 b^i = 2 ^ 6 = 64 
 b^i = 2 ^ 7 = 128 
 b^i = 2 ^ 8 = 256 
 b^i = 2 ^ 9 = 512 
 b^i = 2 ^ 10 = 1024 
 b^i = 2 ^ 11 = 2048 
 b^i = 2 ^ 12 = 4096 
 b^i = 2 ^ 13 = 8192 
 b^i = 2 ^ 14 = 16384 
 b^i = 2 ^ 15 = 32768 
 b^i = 2 ^ 16 = 65536 
 b^i = 2 ^ 17 = 131072 
 b^i = 2 ^ 18 = 262144 
 b^i = 2 ^ 19 = 524288 
 b^i = 2 ^ 20 = 1048576 
 b^i = 2 ^ 21 = 2097152 
 b^i = 2 ^ 22 = 4194304 
 b^i = 2 ^ 23 = 8388608 
 b^i = 2 ^ 24 = 16777216 
 b^i = 2 ^ 25 = 33554432 
 b^i = 2 ^ 26 = 67108864 
 b^i = 2 ^ 27 = 134217728 
 b^i = 2 ^ 28 = 268435456 
 b^i = 2 ^ 29 = 536870912 
 b^i = 2 ^ 30 = 1073741824 
 b^i = 2 ^ 31 = -2147483648 
 b^i = 2 ^ 32 = 0 
 b^i = 2 ^ 33 = 0 
 b^i = 2 ^ 34 = 0 
 b^i = 2 ^ 35 = 0 
 b^i = 2 ^ 36 = 0 
 b^i = 2 ^ 37 = 0 
 b^i = 2 ^ 38 = 0 
 b^i = 2 ^ 39 = 0

Na podstawie wyniku możemy ocenić że zmienna int jest typu 32 bitowego , ponieważ obliczenia do 2^30 są poprawne.

Dla liczb bez znaku przekroczenie zakresu powoduje inny efekt ( modulo ) :

#include <stdio.h>

/* 

 
 Unsigned integers are defined to wrap around.
 "When you work with unsigned types, modular arithmetic (also known as "wrap around" behavior) is taking place."
http://stackoverflow.com/questions/7221409/is-unsigned-integer-subtraction-defined-behavior

*/


int main() {

unsigned int i;
unsigned int old=0; // 
unsigned int new=0; // 
unsigned int p=1000000000; // 
//
unsigned long long int lnew= 0; //
unsigned long long int lold = (unsigned long long int) old; //
unsigned long long int lp = (unsigned long long int) p; //

printf("unsigned long long int \tunsigned  int \n"); // header 

for ( i=0 ; i<20; i++){
     printf("lnew = %12llu \tnew =  %12u",  lnew, new);
     // check overflow
     // http://stackoverflow.com/questions/2633661/how-to-check-integer-overflow-in-c/
     if ( new < old) printf("    unsigned integer overflow = wrap \n");
                     else printf("\n");          
                       
                        
     // unsigned int
     old=new; // save old value for comparison = overflow check 
     new = old + p ; // simple addition ; new value should be greater then old value 
     // unsigned long long int
     lold=lnew;
     lnew=lold+lp;      
}


return 0; 
}

Wynik :

unsigned long long int 	unsigned  int 
lnew =            0 	new =             0
lnew =   1000000000 	new =    1000000000
lnew =   2000000000 	new =    2000000000
lnew =   3000000000 	new =    3000000000
lnew =   4000000000 	new =    4000000000
lnew =   5000000000 	new =     705032704    unsigned integer overflow = wrap 
lnew =   6000000000 	new =    1705032704
lnew =   7000000000 	new =    2705032704
lnew =   8000000000 	new =    3705032704
lnew =   9000000000 	new =     410065408    unsigned integer overflow = wrap 
lnew =  10000000000 	new =    1410065408
lnew =  11000000000 	new =    2410065408
lnew =  12000000000 	new =    3410065408
lnew =  13000000000 	new =     115098112    unsigned integer overflow = wrap 
lnew =  14000000000 	new =    1115098112
lnew =  15000000000 	new =    2115098112
lnew =  16000000000 	new =    3115098112
lnew =  17000000000 	new =    4115098112
lnew =  18000000000 	new =     820130816    unsigned integer overflow = wrap 
lnew =  19000000000 	new =    1820130816

Zapobieganie

sprawdzanie danych :^[30]
- przed wykonaniem działań ^[31]^[32]
  - użycie makr^[33]
- po wykonaniu działań ( może być niebezpieczne dla liczb ze znakiem ponieważ niezdefiniowane zachowanie zagraża bezpieczeństwu komputera )
zwiększenie limitów poprzez :
- zmianę typu ( int , long int, long long int )
- użycie biblioteki o dowolnej precyzji ( np. GMP )

Zapobieganie: wykrywanie możliwego przepełnienia przed wykonaniem działania. Porównaj:

kod z scaler
kod z c-FAQ^[34]

rozmiar

/*
Here is a small C program 
that will print out the size in bytes 
of some basic C types on your machine. 


Paul Gribble | Summer 2012
This work is licensed under a Creative Commons Attribution 4.0 International License
http://gribblelab.org/CBootcamp/3_Basic_Types_Operators_And_Expressions.html

gcc b.c -Wall
./a.out
*/


#include <stdio.h>

int main(int argc, char *argv[]) {
        printf("a char is %ld bytes\n", sizeof(char));
        printf("an int is %ld bytes\n", sizeof(int));
        printf("an float is %ld bytes\n", sizeof(float));
        printf("a double is %ld bytes\n", sizeof(double));
        printf("a short int is %ld bytes\n", sizeof(short int));
        printf("a long int is %ld bytes\n", sizeof(long int));
        printf("a long double is %ld bytes\n", sizeof(long double));
        return 0;
}

a char is 1 bytes
an int is 4 bytes
an float is 4 bytes
a double is 8 bytes
a short int is 2 bytes
a long int is 8 bytes
a long double is 16 bytes

Liczba cyfr

Liczba cyfr w liczbie zmiennoprzecinkowej ^[35]

// http://ubuntuforums.org/showthread.php?t=986212
// http://www.cplusplus.com/reference/cfloat/
// gcc d.c -lm -Wall
// ./a.out

#include <stdio.h>
#include <float.h> 

int main(void)
{
    printf("Float can ensure %d decimal places\n", FLT_DIG);
    printf("Double can ensure %d decimal places\n", DBL_DIG);
    printf("Long double can ensure %d decimal places\n", LDBL_DIG);

    return 0;
}

Wynik :

Float can ensure 6 decimal places
Double can ensure 15 decimal places
Long double can ensure 18 decimal places

Liczby subnormalne

przybliżenia DBL_MIN i liczby subnormalnej

Korzystając z funkcji isnormal zdefiniowanej w pliku math.h możemy samodzielnie poszukać przybliżenia DBL_MIN i liczby subnormalnej.

/* 
isnormal  example 
ISO C99
http://www.cplusplus.com/reference/cmath/isnormal/
http://www.gnu.org/software/libc/manual/html_node/Floating-Point-Classes.html
http://docs.oracle.com/cd/E19957-01/806-3568/ncg_math.html

compile with: 
gcc -std=c99 s.c -lm

run :
./a.out

*/

#include <stdio.h> /* printf */
#include <math.h> /* isnormal */

int TestNumber(double x)
{
  int f; // flag

  f= isnormal(x);
  if (f) 
    printf (" = %g ; number is normal \n",x);
    else printf (" = %g ; number is not normal = denormal  \n",x);

  return f;
}


int main()
{

  double d ;
  double MinNormal; 
  int flag;
  
  d = 1.0 ; // normal
  flag = TestNumber(d);
  do 
  { 
   MinNormal=d;
   d /=2.0; // look for subnormal 
   flag = TestNumber(d);
   
   }
  while (flag);

  printf ("number %f = %g = %e is a approximation of minimal positive double normal \n",MinNormal, MinNormal, MinNormal);
  printf ("number %f = %g = %e is not normal ( subnormal) \n",d, d , d);
   
 return 0;
}

Wynik działania :

number 0.000000 = 2.22507e-308 = 2.225074e-308 is a approximation of minimal positive double normal 
number 0.000000 = 1.11254e-308 = 1.112537e-308 is not normal ( subnormal)

eliminacja liczb subnormalnych

Ten program generuje liczby subnormale:

/*
https://blogs.oracle.com/d/subnormal-numbers
gcc -O0 f.c 
*/
#include <stdio.h>
void main()
{
  double d=1.0;
  while (d>0) {printf("%e\\n",d); d=d/2.0;}
}

wynik:

1.000000e+00
5.000000e-01
2.500000e-01
1.250000e-01
6.250000e-02
3.125000e-02
...
3.162020e-322
1.581010e-322
7.905050e-323
3.952525e-323
1.976263e-323
9.881313e-324
4.940656e-324

Jeśli jednak skompilujemy go z opcję:

  gcc -O0 -ffast-math f.c

to otrzymamy:

...
3.560118e-307
1.780059e-307
8.900295e-308
4.450148e-308
2.225074e-308

Liczby subnormalne są zaokrąglane do zera.

Jaki wpływ na obliczenia mają liczby subnormalne?

wydłużają czas obliczeń^[36]

część ułamkowa

Za pomocą:^[37]

funkcji modf
konwersji int

double frac = r - (int)r ;

Błędy w obliczeniach numerycznych

Typy:^[38]^[39]

wg etapu operacji :^[40]
- blędne dane wejściowe : niezgodne z oczekiwanym typem
- dane wejściowe powodują błąd rezultatu
wg rodzaju operacji
- zmiennoprzecinkowych^[41] ^[42]^[43]^[44]
  - porównywanie ^[45]
  - odejmowanie podobnych liczb^[46]^[47]
  - błędy zaokrąglania ^[48]
  - przekroczenie limitu ( underflow, overflow)^[49]^[50]
  - mnożenie

Efekt:

zachowanie niezdefiniowane : ( ang. undefined behavior = UB)^[51] ^[52]

Dzielenie przez zero

Dzielenie przez zero^[53]

Kiedy dzielnik ma wartość zero w wyrażeniu (zwykle przez pomyłkę) to następuje awaria programu (nieprawidłowe zakończenie an. crash). Nieprawidłowe zakończenie może być poważnym problemem w przypadku oprogramowania o krytycznym znaczeniu dla życia.

Zapobiegać temu można przez:

kontrole if-else
obsługę wyjątków

Mnożenie

obliczanie silni

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
/* 
https://math.stackexchange.com/questions/2453939/is-this-characteristic-of-tent-map-usually-observed
*/

/* ------------ constans  ---------------------------- */
double m = 2.0; /* parameter of tent map */
double a = 1.0; /* upper bound for randum number generator */
int iMax = 100;
/* ------------------- functions --------------------------- */

/* 
tent map 
https://en.wikipedia.org/wiki/Tent_map
*/
double f(double x0, double m){

	double x1;
	if (x0 < 0.5) 
		x1 = m*x0;
		else x1 = m*(1.0 - x0);
	return x1;

}



/* random double from 0.0 to a 
https://stackoverflow.com/questions/13408990/how-to-generate-random-float-number-in-c


*/
double GiveRandom(double a){
	srand((unsigned int)time(NULL));
	
	return  (((double)rand()/(double)(RAND_MAX)) * a);

}

int main(void){

	int i = 0;
	double x = GiveRandom(a); /* x0 = random */
	
	for (i = 0; i<iMax; i++){
	
		printf("i = %3d \t x = %.16f\n",i, x);
		x = f(x,m); /* iteration of the tent map */
	}

	return 0;
}

Kompilacja i uruchomienie:

  gcc t.c -Wall
  ./a.out

Wynik:

i =   0 	 x = 0.1720333817284710
i =   1 	 x = 0.3440667634569419
i =   2 	 x = 0.6881335269138839
i =   3 	 x = 0.6237329461722323
i =   4 	 x = 0.7525341076555354
i =   5 	 x = 0.4949317846889292
i =   6 	 x = 0.9898635693778584
i =   7 	 x = 0.0202728612442833
i =   8 	 x = 0.0405457224885666
i =   9 	 x = 0.0810914449771332
i =  10 	 x = 0.1621828899542663
i =  11 	 x = 0.3243657799085327
i =  12 	 x = 0.6487315598170653
i =  13 	 x = 0.7025368803658694
i =  14 	 x = 0.5949262392682613
i =  15 	 x = 0.8101475214634775
i =  16 	 x = 0.3797049570730451
i =  17 	 x = 0.7594099141460902
i =  18 	 x = 0.4811801717078197
i =  19 	 x = 0.9623603434156394
i =  20 	 x = 0.0752793131687213
i =  21 	 x = 0.1505586263374425
i =  22 	 x = 0.3011172526748851
i =  23 	 x = 0.6022345053497702
i =  24 	 x = 0.7955309893004596
i =  25 	 x = 0.4089380213990808
i =  26 	 x = 0.8178760427981615
i =  27 	 x = 0.3642479144036770
i =  28 	 x = 0.7284958288073540
i =  29 	 x = 0.5430083423852921
i =  30 	 x = 0.9139833152294159
i =  31 	 x = 0.1720333695411682
i =  32 	 x = 0.3440667390823364
i =  33 	 x = 0.6881334781646729
i =  34 	 x = 0.6237330436706543
i =  35 	 x = 0.7525339126586914
i =  36 	 x = 0.4949321746826172
i =  37 	 x = 0.9898643493652344
i =  38 	 x = 0.0202713012695312
i =  39 	 x = 0.0405426025390625
i =  40 	 x = 0.0810852050781250
i =  41 	 x = 0.1621704101562500
i =  42 	 x = 0.3243408203125000
i =  43 	 x = 0.6486816406250000
i =  44 	 x = 0.7026367187500000
i =  45 	 x = 0.5947265625000000
i =  46 	 x = 0.8105468750000000
i =  47 	 x = 0.3789062500000000
i =  48 	 x = 0.7578125000000000
i =  49 	 x = 0.4843750000000000
i =  50 	 x = 0.9687500000000000
i =  51 	 x = 0.0625000000000000
i =  52 	 x = 0.1250000000000000
i =  53 	 x = 0.2500000000000000
i =  54 	 x = 0.5000000000000000
i =  55 	 x = 1.0000000000000000
i =  56 	 x = 0.0000000000000000
i =  57 	 x = 0.0000000000000000
i =  58 	 x = 0.0000000000000000
i =  59 	 x = 0.0000000000000000
i =  60 	 x = 0.0000000000000000
i =  61 	 x = 0.0000000000000000
i =  62 	 x = 0.0000000000000000
i =  63 	 x = 0.0000000000000000
i =  64 	 x = 0.0000000000000000
i =  65 	 x = 0.0000000000000000
i =  66 	 x = 0.0000000000000000
i =  67 	 x = 0.0000000000000000
i =  68 	 x = 0.0000000000000000
i =  69 	 x = 0.0000000000000000
i =  70 	 x = 0.0000000000000000
i =  71 	 x = 0.0000000000000000
i =  72 	 x = 0.0000000000000000
i =  73 	 x = 0.0000000000000000
i =  74 	 x = 0.0000000000000000
i =  75 	 x = 0.0000000000000000
i =  76 	 x = 0.0000000000000000
i =  77 	 x = 0.0000000000000000
i =  78 	 x = 0.0000000000000000
i =  79 	 x = 0.0000000000000000
i =  80 	 x = 0.0000000000000000
i =  81 	 x = 0.0000000000000000
i =  82 	 x = 0.0000000000000000
i =  83 	 x = 0.0000000000000000
i =  84 	 x = 0.0000000000000000
i =  85 	 x = 0.0000000000000000
i =  86 	 x = 0.0000000000000000
i =  87 	 x = 0.0000000000000000
i =  88 	 x = 0.0000000000000000
i =  89 	 x = 0.0000000000000000
i =  90 	 x = 0.0000000000000000
i =  91 	 x = 0.0000000000000000
i =  92 	 x = 0.0000000000000000
i =  93 	 x = 0.0000000000000000
i =  94 	 x = 0.0000000000000000
i =  95 	 x = 0.0000000000000000
i =  96 	 x = 0.0000000000000000
i =  97 	 x = 0.0000000000000000
i =  98 	 x = 0.0000000000000000
i =  99 	 x = 0.0000000000000000

Porównywanie

Sprawdźmy czy liczba x jest równa zero :

if (x==0.0)

Czy takie porównanie jest bezpieczne dla liczb zmiennoprzecinkowych ?^[54]

//  gcc c1.c -Wall -lm
#include <math.h> /* isnormal */
#include <float.h>//DBL_MIN
#include <stdio.h>

int main ()
{
  double x = 1.0;
 

  
 
  int i;
  for ( i=0; i < 334; i++)
    {
      x/=10;
      printf ("i = %3d ; x= %.16lf = %e so  ", i, x,x);
      //
      if (x<DBL_MIN) printf ("x < DBL_MIN and ");
          else printf ("x > DBL_MIN  and  ");
      //
      if (isnormal(x)) printf ("x is normal and ");
          else printf ("x is  subnormal and  ");
      //
      if (x==0.0) printf ("equal to 0.0\n");
      else printf ("not equal to 0.0\n");
    }

return 0;
}

Wynik :

i =   0 ; x= 0.1000000000000000 = 1.000000e-01 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   1 ; x= 0.0100000000000000 = 1.000000e-02 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   2 ; x= 0.0010000000000000 = 1.000000e-03 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   3 ; x= 0.0001000000000000 = 1.000000e-04 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   4 ; x= 0.0000100000000000 = 1.000000e-05 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   5 ; x= 0.0000010000000000 = 1.000000e-06 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   6 ; x= 0.0000001000000000 = 1.000000e-07 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   7 ; x= 0.0000000100000000 = 1.000000e-08 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   8 ; x= 0.0000000010000000 = 1.000000e-09 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =   9 ; x= 0.0000000001000000 = 1.000000e-10 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  10 ; x= 0.0000000000100000 = 1.000000e-11 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  11 ; x= 0.0000000000010000 = 1.000000e-12 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  12 ; x= 0.0000000000001000 = 1.000000e-13 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  13 ; x= 0.0000000000000100 = 1.000000e-14 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  14 ; x= 0.0000000000000010 = 1.000000e-15 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  15 ; x= 0.0000000000000001 = 1.000000e-16 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  16 ; x= 0.0000000000000000 = 1.000000e-17 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  17 ; x= 0.0000000000000000 = 1.000000e-18 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  18 ; x= 0.0000000000000000 = 1.000000e-19 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  19 ; x= 0.0000000000000000 = 1.000000e-20 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i =  20 ; x= 0.0000000000000000 = 1.000000e-21 so  x > DBL_MIN  and  x is normal and not equal to 0.0
...
i = 290 ; x= 0.0000000000000000 = 1.000000e-291 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 291 ; x= 0.0000000000000000 = 1.000000e-292 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 292 ; x= 0.0000000000000000 = 1.000000e-293 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 293 ; x= 0.0000000000000000 = 1.000000e-294 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 294 ; x= 0.0000000000000000 = 1.000000e-295 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 295 ; x= 0.0000000000000000 = 1.000000e-296 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 296 ; x= 0.0000000000000000 = 1.000000e-297 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 297 ; x= 0.0000000000000000 = 1.000000e-298 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 298 ; x= 0.0000000000000000 = 1.000000e-299 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 299 ; x= 0.0000000000000000 = 1.000000e-300 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 300 ; x= 0.0000000000000000 = 1.000000e-301 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 301 ; x= 0.0000000000000000 = 1.000000e-302 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 302 ; x= 0.0000000000000000 = 1.000000e-303 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 303 ; x= 0.0000000000000000 = 1.000000e-304 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 304 ; x= 0.0000000000000000 = 1.000000e-305 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 305 ; x= 0.0000000000000000 = 1.000000e-306 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 306 ; x= 0.0000000000000000 = 1.000000e-307 so  x > DBL_MIN  and  x is normal and not equal to 0.0
i = 307 ; x= 0.0000000000000000 = 1.000000e-308 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 308 ; x= 0.0000000000000000 = 1.000000e-309 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 309 ; x= 0.0000000000000000 = 1.000000e-310 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 310 ; x= 0.0000000000000000 = 1.000000e-311 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 311 ; x= 0.0000000000000000 = 1.000000e-312 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 312 ; x= 0.0000000000000000 = 1.000000e-313 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 313 ; x= 0.0000000000000000 = 1.000000e-314 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 314 ; x= 0.0000000000000000 = 1.000000e-315 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 315 ; x= 0.0000000000000000 = 1.000000e-316 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 316 ; x= 0.0000000000000000 = 9.999997e-318 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 317 ; x= 0.0000000000000000 = 9.999987e-319 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 318 ; x= 0.0000000000000000 = 9.999889e-320 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 319 ; x= 0.0000000000000000 = 9.999889e-321 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 320 ; x= 0.0000000000000000 = 9.980126e-322 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 321 ; x= 0.0000000000000000 = 9.881313e-323 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 322 ; x= 0.0000000000000000 = 9.881313e-324 so  x < DBL_MIN and x is  subnormal and  not equal to 0.0
i = 323 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 324 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 325 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 326 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 327 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 328 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 329 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 330 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 331 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 332 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0
i = 333 ; x= 0.0000000000000000 = 0.000000e+00 so  x < DBL_MIN and x is  subnormal and  equal to 0.0

Jak powinno się porównywać liczby zmienno przecinkowe ?^[55]

wartość bezwględna róznicy : if( abs(a-b) < epsilon) // wrong - don't do this^[56]
if( abs((a-b)/b) < epsilon ) // still not right!
wartości graniczne
stałe ^[57]

 if (fpclassify(x) == FP_ZERO )

lub

 if (x == FP_ZERO)

Wartości służące do testo wania porównań :

wg Michael Borgwardt^[58]

Sumowanie

Na ile poważny jest to problem? Spróbujmy przyjrzeć się działaniu, polegającym na 1000-krotnym dodawaniu do liczby wartości 1/3. Oto kod:

#include <stdio.h>

int main ()
{
  float a = 0;
  int i = 0;
  for (;i<1000;i++)
  {
    a += 1.0/3.0;
  }
  printf ("%f\n", a);
}

Z matematyki wynika, że 1000*(1/3) = 333.(3), podczas gdy komputer wypisze wynik, nieco różniący się od oczekiwanego (w moim przypadku):

333.334106

Błąd pojawił się na cyfrze części tysięcznej liczby. Nie jest to może poważny błąd, jednak zastanówmy się, czy ten błąd nie będzie się powiększał. Zamieniamy w kodzie ilość iteracji z 1000 na 100 000. Tym razem mój komputer wskazał już nieco inny wynik:

33356.554688

Błąd przesunął się na cyfrę dziesiątek w liczbie. Tak więc nie należy do końca polegać na prezentacji liczb zmiennoprzecinkowych w komputerze.

Utrata precyzji

Utrata precyzji, utrata cyfr znaczących ( ang. Loss of significance, catastrophic cancellation of significant digits)

sumowanie dużej liczby z małą
odejmowanie prawie równych liczb^[59]

Biblioteki matematyczne

float.h

opis

math.h

Aby móc korzystać z wszystkich dobrodziejstw funkcji matematycznych musimy na początku dołączyć plik math.h:

 #include <math.h>

a w procesie kompilacji (dotyczy kompilatora GCC) musimy dodać flagę "-lm" po nazwie pliku wynikowego^[60], czyli na końcu linii :^[61]

Flaga lm jest zależna od środowiska. Na przykład, w systemie Windows, nie jest to wymagana, ale jest to wymagana w systemach opartych na systemie UNIX.

gcc plik.c -o plik -lm

Funkcje matematyczne, które znajdują się w bibliotece standardowej ( plik libm.a ) możesz znaleźć tutaj. Przy korzystaniu z nich musisz wziąć pod uwagę m.in. to, że biblioteka matematyczna prowadzi kalkulację w oparciu o radiany a nie stopnie.

Wersje

Stałe matematyczne: pi, e ...

W pliku math.h zdefiniowane są pewne stałe, które mogą być przydatne do obliczeń. Są to m.in.:

M_E - podstawa logarytmu naturalnego (e, liczba Eulera)
M_LOG2E - logarytm o podstawie 2 z liczby e
M_LOG10E - logarytm o podstawie 10 z liczby e
M_LN2 - logarytm naturalny z liczby 2
M_LN10 - logarytm naturalny z liczby 10
M_PI - liczba π
M_PI_2 - liczba π/2
M_PI_4 - liczba π/4
M_1_PI - liczba 1/π
M_2_PI - liczba 2/π

Możemy to sprawdzić:

grep -i pi /usr/include/math.h

i otrzymamy:

# define M_PI        3.14159265358979323846    /* pi */
# define M_PI_2        1.57079632679489661923    /* pi/2 */
# define M_PI_4        0.78539816339744830962    /* pi/4 */
# define M_1_PI        0.31830988618379067154    /* 1/pi */
# define M_2_PI        0.63661977236758134308    /* 2/pi */
# define M_2_SQRTPI    1.12837916709551257390    /* 2/sqrt(pi) */
# define M_PIl        3.1415926535897932384626433832795029L  /* pi */
# define M_PI_2l    1.5707963267948966192313216916397514L  /* pi/2 */
# define M_PI_4l    0.7853981633974483096156608458198757L  /* pi/4 */
# define M_1_PIl    0.3183098861837906715377675267450287L  /* 1/pi */
# define M_2_PIl    0.6366197723675813430755350534900574L  /* 2/pi */
# define M_2_SQRTPIl    1.1283791670955125738961589031215452L  /* 2/sqrt(pi) */
/* When compiling in strict ISO C compatible mode we must not use the

Liczby zespolone

teoria ( wikibooks)

Typy

_Complex
complex.h
biblioteki:
- mpc
- arb ^[62]

różnica pomiędzy _complex a complex

_Complex jest słowem kluczowym, możemy go używać bez dyrektywy include dołączającej plik nagłówkowy complex.h
complex jest makrem z pliku nagłówkowego complex.h

float _Complex     
double _Complex                  
long double _Complex

 #include <complex.h>
 float complex
double complex
long double complex

complex.h

Operacje na liczbach zespolonych są częścią uaktualnionego standardu języka C o nazwie C99, który jest obsługiwany jedynie przez część kompilatorów

Dotychczas korzystaliśmy tylko z liczb rzeczywistych, lecz najnowsze standardy języka C umożliwiają korzystanie także z innych liczb - np. z liczb zespolonych.

Aby móc korzystać z liczb zespolonych w naszym programie należy w nagłówku programu umieścić następującą linijkę:

 #include <complex.h>

która powoduje dołączenie standardowej biblioteki obsługującej liczny zespolenie

Wiemy, że liczba zespolona zdeklarowana jest następująco:

z = a+b*i,

gdzie a, b są liczbami rzeczywistymi, a i jest jednostką urojoną

i*i = (-1).

W pliku complex.h liczba i zdefiniowana jest jako I. Zatem wypróbujmy możliwości liczb zespolonych:

 #include <math.h>
 #include <complex.h>
 #include <stdio.h>
 
 int main ()
 {
   float _Complex z = 4+2.5*I;
   printf ("Liczba z: %f+%fi\n", creal(z), cimag (z));
   return 0;
 }

następnie kompilujemy nasz program:

gcc plik1.c -o plik1 -lm

Po wykonaniu naszego programu powinniśmy otrzymać:

Liczba z: 4.00+2.50i

W programie zamieszczonym powyżej użyliśmy dwóch funkcji - creal i cimag.

creal - zwraca część rzeczywistą liczby zespolonej
cimag - zwraca część urojoną liczby zespolonej

Więcej:

// https://stackoverflow.com/questions/6418807/how-to-work-with-complex-numbers-in-c
// program by user870774
#include <stdio.h>      /* Standard Library of Input and Output */
#include <complex.h>    /* Standard Library of Complex Numbers */

int main() {

    double complex z1 = 1.0 + 3.0 * I;
    double complex z2 = 1.0 - 4.0 * I;

    printf("Working with complex numbers:\n\v");

    printf("Starting values: Z1 = %.2f + %.2fi\tZ2 = %.2f %+.2fi\n", creal(z1), cimag(z1), creal(z2), cimag(z2));

    double complex sum = z1 + z2;
    printf("The sum: Z1 + Z2 = %.2f %+.2fi\n", creal(sum), cimag(sum));

    double complex difference = z1 - z2;
    printf("The difference: Z1 - Z2 = %.2f %+.2fi\n", creal(difference), cimag(difference));

    double complex product = z1 * z2;
    printf("The product: Z1 x Z2 = %.2f %+.2fi\n", creal(product), cimag(product));

    double complex quotient = z1 / z2;
    printf("The quotient: Z1 / Z2 = %.2f %+.2fi\n", creal(quotient), cimag(quotient));

    double complex conjugate = conj(z1);
    printf("The conjugate of Z1 = %.2f %+.2fi\n", creal(conjugate), cimag(conjugate));

    return 0;
}

MPC

Opis MPC

Dodatkowe

Zobacz również

Źródła

« Więcej o kompilowaniu

Spis treści

Powszechne praktyki »

[1] The Rational Number Class in C by Peter Burden

[2] Rational Arithmetic by R. Sedgewick

[3] Where is the itoa function in Linux?

[4] - Binary-constants

[5] Code Replay : C represent int in base 2

[6] what-every-computer-programmer-should by Josh Haberman

[7] Metody numeryczne - autorzy : Piotr Krzyżanowski i Leszek Plaskota — Uniwersytet Warszawski, Wydział Matematyki, Informatyki i Mechaniki

[8] How to get memory representation-double

[9] IEEE-754 Floating Point Converter

[10] umpfp: A Tool to Inspect Floating-Point Numbers by Joshua Haberman

[11] subethasoftware : converting-two-8-bit-values-to-one-16-bit-value-in-c

[12] subethasoftware: splitting-a-16-bit-value-to-two-8-bit-values-in-c

[13] Math to Code by cameron smith

[14] veloping-mathematical-software-in-c by Fredrik Johansson

[15] rinting-algebraic-numbers by Fredrik Johansson

[16] umeryczne - Wydziału MIM UW

[17] Cpp Core Guidelines : arithmetic

[18] ól obliczeń numerycznych -Piotr Krzyżanowski

[19] Every Computer Scientist Should Know About Floating-Point Arithmetic, by David Goldberg

[20] Two disasters caused by computer arithmetic errors

[21] Praktyczne wyznaczanie precyzji arytmetyki - autorzy : Piotr Krzyżanowski i Leszek Plaskota — Uniwersytet Warszawski, Wydział Matematyki, Informatyki i Mechaniki

[22] stackoverflow question : why-does-gcc-report-a-floating-point-exception-when-i-execute-1-0

[23] stackoverflow question: catch-floating-point-exceptions-using-a-compiler-option-with-c

[24] Point Representation - Basics from : geeksforgeeks.org

[25] Tutorial on Data Representation by Chua Hock-Chuan

[26] Przekroczenie zakresu liczb całkowitych w wikipedii

[27] [http://www.fefe.de/intof.html%7CCatching Integer Overflows in C ]

[28] Guide to Undefined Behavior in C and C++, Part 1 by John Regehr

[29] stackoverflow question: function-abs-returning-negative-number-in-c

[30] Testing for Integer Overflow in C and C++ by Josh Haberman

[31] .lang.c FAQ list · Question 20.6b : How can I ensure that integer arithmetic doesn't overflow ?

[32] GCC : Built-in Functions to Perform Arithmetic with Overflow Checking

[33] ulib : Checking-Integer-Overflow

[34] -faq: three functions for ``careful addition, subtraction, and multiplication

[35] Stackoverflow : Counting digits in a float

[36] O N SUBNORMAL FLOATING POINT AND ABNORMAL TIMING by Marc Andrysco, David Kohlbrenner, Keaton Mowery, Ranjit Jhala, Sorin Lerner, and Hovav Shacham

[37] stackoverflow question : extract-decimal-part-from-a-floating-point-number-in-c

[38] INTRODUCTION TO NUMERICAL ANALYSIS WITH C PROGRAMS by Attila Mate

[39] Basic Issues in Floating Point Arithmetic and Error Analysis by Jim Demmel

[40] securecoding.cert.org : Prevent+or+detect+domain+and+range+errors+in+math+functions

[41] Floating point inaccuracy examples

[42] Catastrophic Cancellation: The Pitfalls of Floating Point Arithmetic by Graham Markall

[43] th.stackexchange question: is-this-characteristic-of-tent-map-usually-observed

[44] BŁĘDY PRZETWARZANIA NUMERYCZNEGO, Maciej Patan, UW Zielonogórski

[45] rous is it to compare floating point values?

[46] t possible to get 0 by subtracting two unequal floating point numbers?

[47] How to solve quadratic equations numerically by FLORIAN DANG

[48] Double Rounding Errors in Floating-Point Conversions By Rick Regan (Published August 11th, 2010)

[49] stackoverflow question when-does-underflow-occur

[50] stackoverflow question: how-to-define-underflow-for-an-implementationieee754-which-support-subnormal-n

[51] Undefined_behavior w ang. wikipedii

[52] undefined-behavior-in-c-and-cplusplus-programs by Nayuki

[53] tstack :infinity-in-c by Abdul Mateen Oct-26, 2022

[54] stackoverflow question: c-floating-point-zero-comparison

[55] Comparing Floating-Point Numbers Is Tricky by Matt Kline

[56] | floating-point-gui.de : comparison/

[57] stackoverflow question : float-double-equality-with-exact-zero

[58] Michael Borgwardt : Nearly Equals Test in java

[59] umerically-stable-law-of-cosines by Nayuki

[60] Stackoverflow : Undefined reference to `pow' and `floor'

[61] I'm using math.h and the library link -lm, but “undefined reference to `pow'” still happening

[62] rb library

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]

[11]

[12]

[13]

[14]

[15]

[16]

[17]

[18]

[19]

[20]

[21]

[22]

[23]

[24]

[25]

[26]

[27]

[28]

[29]

[30]

[31]

[32]

[33]

[34]

[35]

[36]

[37]

[38]

[39]

[40]

[41]

[42]

[43]

[44]

[45]

[46]

[47]

[48]

[49]

[50]

[51]

[52]

[53]

[54]

[55]

[56]

[57]

[58]

[59]

[60]

[61]

[62]