Převod iracionálního čísla na zlomek - VYŘEŠENO

LADER · 9. 11. 2014 10:30:52

Zdravím.
Potřebuji převést iracionální číslo na poměr dvou celých malých čísel, které co nejpřesněji vyjadřují původní číslo.

Příklad:
1.713933142811 převede na 12/7
nebo:
sqr(5-sqr(2)) se převede na 89/47

Nejlépe tak, aby čitatel byl v jedné buňce a jmenovatel ve druhé.

Ještě poznámka, nechci aby funkce hledala zlomek podle požadované přesnosti, protože pak zlomek sice může vyhovovat dané přesnosti, ale nemusí vyhovovat na jednoduchost.

Příklad:
číslo Pí s přesností lepší než 10^-5 vyjde 3927/1250,
ale je lepší možnost 355/113, které má podstatně vyšší přesnost.

Takže bych chtěl, aby ta funkce vypisovala pouze ta "výhodná" řešení. Přesnost se dá dopočítat dodatečně.

Neznáte nějakou hotovou jednoduchou funkci (případně funkční makro)?

Editoval LADER (9. 11. 2014 16:26:08)

neutr · 9. 11. 2014 12:43:03

Takováhle úloha patří do teorie čísel. Nicméně lze si pomoci. Zaokrouhlíme iracionální číslo na určitý počet desetinných míst. Vztah mezi dělencem a dělitelem je dán samozřejmě jako podíl, ale také je to můžné vyjádřit jako rovnici pomocí MODULO a INTEGER Totiž :
INTEGER
(INT) = Vzorcem =INT(Dělenec;Dělitel) - Basic : INT(Dělenec/dělitel)
MODULO
(MOD) = Vzorcem =MOD(Dělenec;Dělitel) - Basic : Dělenec MOD dělitel

Zlomek obecně potom ve tvaru desetinného čísla odpovídá následujícímu vztahu například :
61,2577724888
- Celočíselná část (červená) odpovídá jak výsledek INTEGER
- Neceločíselná část (modrá) je v relaci s MODULO (Modulo je celočíselný zbytek dělenec - Integer)
Zde se může postup rozvětvit (více možných řešení). Základem je iterace možných dělitelů (D) M(D)[2...(d-1)]..(d)=dělenec.
Takže iterace pomocí vzorce by také přicházela do úvahy spolu s řešitelem. Je otázkou nastavení systémových iterací a iterací potřebných v rozsahu [2...(d-1)]. Ovšem pro relativně malá čísla by to mohlo vyhovovat.

Pokud by to bylo opravdu málo možností (malá celočíselná část), musí se vlastní podíl přiměřeně zvětšovat. Jde li tedy o podíl < 1, je nutno pracovat s jeho převrácenou hodnnotou.

Zde by mohly výrazně pomoci maticové vzorce. Což by byla asi jiná alternativa postupu, nežli s těmi systémovými iteracemi.

Samozřejmě zde přichází nejspíš do úvahy makro, respektive psaná funkce, která se spouští ze sešitu jako vzorec. Tady bych to viděl jako práci s array.

Jen pro orientaci : číslo za zlomkovou čárkou (x) je v relaci s Modulo takto :
"x" = MODULO / dělenec => dělenec = MODULO/(x) - respektive pro lepší představu (celočíselný zbytek po dělení/dělitel).
Potom hledaný celočíselný dělenec (d) a dělitel (D) :
(d) = INT*(i)
(D) = MOD/(i)
Odhad (D) = cca 1/(x) - tedy převrácené hodnotě neceločíselné části. Z toho dostaneme startovní hodnoty pro iterace podílů. jestliže 1/(x) = (y), potom :
Celočíselná část desetinného čísla (A) * (y) = cca hledanému čitateli (Č).
Převrácená hodnota (1/x) = cca jmenovatel (J).

Dle mne nejlepší postup : Vycházáme z toho, že (Č) + (J) = číslo z oboru čísel N. Rozklad tohoto čísla podle Partition Numerorum[PN] určí posloupnost dělenců a dělitelů.
Přibližovat se můžeme jednoduše pomocí inkrementace , nebo dekrementace (Č + J)+/-1. Vlastní hledání pak musí vycházet z toho zda a jak velká je celočíselná část desetinného čísla.
Například třída PN(5,2) vypadá takto :
5 + 0
4 + 1
3 + 2
2 + 3
1 + 4
0 + 5
Je logické, že desetinné číslo < 1 budeme hledat mexi 2/3 a 1/4. Naopak číslo > 1 mezi 4/1 a 3/2. Pakliže je to nevyhovující provedeme iterace na PN(6,2).

Algoritmus je jednoduchý a iterace zpracujeme pomocí array - což umíte. Jen bude asi potřeba někdy zpracovat výpočtovou iteraci + jednu dekremenovanou a také inkrementovanou.
Stejně bychom postupovali při určování přesnosti. Takže bych to viděl na funkci spouštěnou ze sešitu, které dáte dva parametry - vlastní desetinné číslo a údaj vyjadřující přesnost například v procentech ap.

Je zajímavé k čemu by to bylo potřeba - že by třeba k určování dědických podílů nebo jiných majetkových? Ale samozřejmě mimo justiční praxi se příkladů najde přehršel.
Pokud by bylo potřeba vyšších tříd PN, bez pomoci se asi neobejdete. PN(n,2) - druhá třída je snadná "A" jako dělenec a "B" jako dělitel jsou dány "A"[0..n], "B" = n-"A". Respektive hodnoty nulových dělenců a dělitelů jsou nesmyslné a stejně tak dělenců a dělitelů jedničkových. Takže pro
čísla < 1 je to "A"[2..(1/2n)] a čísla > 1 je to "A"[(1/2n)...(n-1)].

lp. · 9. 11. 2014 13:08:16

Něco jsem našel

Function farey(x As Double, N As Long) As string
' podle http://www.johndcook.com/blog/2010/10/20/best-rational-approximation/
Dim a As Long, b As Long, c As Long, d As Long
Dim xxx As Long
Dim mediant As Double

    xxx = Fix(x)
    x = x - xxx
    
    a = 0
    b = 1
    c = 1
    d = 1
    
    While ((b <= N) And (d <= N))
        mediant = CDbl(a + c) / (b + d)
        If x = mediant Then
            If b + d <= N Then
              farey = xxx & " " & (a + c) & "/" & (b + d)
            ElseIf d > b Then
              farey = xxx & " " & c & "/" & d
            Else
              farey = xxx & " " & a & "/" & b
            End If
            Exit Function
        ElseIf x > mediant Then
            a = a + c
            b = b + d
        Else
            c = a + c
            d = b + d
        End If
 
    Wend
    
    If (b > N) Then
      farey =  xxx & " " & c & "/" & d
    Else
      farey = xxx & " " & a & "/" & b
    End If
End Function

(Varianta s navracení čísle a jejich převodem na zlomek formátem zlomky selhávalo kvůli nepřesnostem při použití tohoto formátu.)

LADER · 9. 11. 2014 13:27:10

Díky za odpověď,
Podle Partition Numerorum[PN] jsem něco napsal, sice to funguje, ale není to ono:

Function Najdi(X As double, E As double) As Variant
' X = zadané kladné iracionální číslo, které se převede na zlomek
' E = požadovaná přesnost (např.: 1E-5)
	Dim i As Integer, N As Double, D As Double, V(0,2) As Double 
	i = 1000  	' Ochrana proti zacyklení
	N = 1		' Čitatel	
	D = 1		' Jmenovatel
	While (i > 0%) And abs(X - N/D) > E
		If X - N/D > 0 Then
			N = N + 1
		Else
			D = D + 1
		EndIf
		i = i - 1%
	Wend
	V(0,0) = N/D
	V(0,1) = N
	V(0,2) = D
	Najdi = V
End Function

Je to pomalé, inkrementuje to pouze o jedničku ...

LADER · 9. 11. 2014 13:46:48

lp. napsal(a)

Něco jsem našel ...

Jo díky, to vypadá velmi dobře. Trochu jsem to přepsal (jenom výstup, ostatní jsem zachoval). Teď to budu testovat.

Function Farey(x As Double, N As Long) As Variant
' Podle http://www.johndcook.com/blog/2010/10/20/best-rational-approximation/
' Upraveno
Dim a As Long, b As Long, c As Long, d As Long
Dim xxx As Long, mediant As Double, V(0,2) As Double
    xxx = Fix(x)
    x = x - xxx    
    a = 0
    b = 1
    c = 1
    d = 1    
    While ((b <= N) And (d <= N))
        mediant = CDbl(a + c) / (b + d)
        If x = mediant Then
            If b + d <= N Then            	
				V(0,1) = (a + c)+xxx*(b + d)
				V(0,2) = (b + d)
				V(0,0) = V(0,1)/V(0,2)
				Farey = V
            ElseIf d > b Then         
				V(0,1) = (c)+xxx*(d)
				V(0,2) = (d)
				V(0,0) = V(0,1)/V(0,2)
				Farey = V
            Else
				V(0,1) = (a)+xxx*(b)
				V(0,2) = (b)
				V(0,0) = V(0,1)/V(0,2)
				Farey = V
            End If
            Exit Function
        ElseIf x > mediant Then
            a = a + c
            b = b + d
        Else
            c = a + c
            d = b + d
        End If 
    Wend
    
    If (b > N) Then
		V(0,1) = (c) + xxx*(d)
		V(0,2) = (d)
		V(0,0) = V(0,1)/V(0,2)
		Farey = V
    Else
		V(0,1) = (a) + xxx*(b)
		V(0,2) = (b)
		V(0,0) = V(0,1)/V(0,2)
		Farey = V
    End If
End Function

Trošku se to rozhodilo (asi jiný tabulátory).

LADER · 9. 11. 2014 16:25:30

Zdravím,
Nová funkce Farey vypadá dobře skutečně vybírá ty nejlepší poměry, ale sem tam najde i méně výhodný poměr (ale pořád ještě dobrý). je to asi tím, že se zadává požadovaná přesnost.
Proto jsem podle Continued fraction vytvořil novou funkci, která tímto problémem netrpí. Zadává se tam počet průchodů a ta každým cyklem najde ten následující nejlepší poměr. Přesnost jde dopočítat dodatečně.
Kód je trochu delší, snad tam nejsou závažné chyby:

Function nMax(X as Double, Y As Double) as Double
	If X>Y Then
		nMax = X
	Else 
		nMax = Y
	End If
End Function

Function nMin(X as Double, Y As Double) as Double
	If X<Y Then
		nMin = X
	Else 
		nMin = Y
	End If
End Function

Function Fraction(X as Double, N As Integer) as Variant
	Dim A(20%) As Long, i As Integer, k As Integer, C As Long, F As Integer
	Dim S As Integer, X1 As Double, Z As Double, D As Double 
	Dim Zn As double, Zd As Double, V(0%,2%) As Double
	Dim Kn(20%) As Double, Kd(20%) As Double
	
	k = 1%
	S = sgn(X) 	' Znaménko
	X = abs(X)	' Absolutní hodnota
	N = nMin(nMax(1%, N), 16%)	' Nastavení minimálního a maximálního počtu opakování
	
	If S <> 0% Then 	' Hodnota nesmí být nulová
		If X < 1 Then 	' Hodnota je menší než 1
			F = 1%
			X = 1# / X
		Else
			F = 0%
		End If
		
		For i = 1% To N	' Rozložení na řetězový zlomek -> výsledek je v poli A			
			C = int(X + 0.5)
			A(i) = C
			k = i
			D = X - C			
			Select Case i	' Vypočítáme čitatele Zn a jmenovatele zlomku Zd (z řetězového zlomku)
			Case 1%
				Kn(1%) = A(1%)
				Kd(1%) = 1#			
			Case 2%
				Kn(2%) = A(1%)*A(2%)+1#
				Kd(2%) = A(2%)				
			Case 3% To 17%
				Kn(i) = A(i)*Kn(i-1%)+Kn(i-2%)
				Kd(i) = A(i)*Kd(i-1%)+Kd(i-2%)					
			End Select
			If abs(D) < 1e-9 Then Exit For		
			X = 1# / D		
		Next i
		Zn = Kn(k)			
		Zd = Kd(k)
						
		If Zd < 0# Then 	' Pokud vyjde záporný jmenovatel, otočíme znaménka
			Zd = -Zd
			Zn = -Zn
		End If
		
		If F = 1% Then	' Pokud bylo převáděné číslo větší jak jedna, vyměníme hodnoty
			D  = Zn
			Zn = Zd
			Zd = D	
		End If 		
	Else	' Pokud byla zadána nula je výsledek nula
		Z = 0#
		Zn = 0#
		Zd = 1#
	End If
	 	
	V(0%,0%) = S * Zn/Zd	' Vypočtený zlomek
	V(0%,1%) = S * Zn		' Numerator
	V(0%,2%) = Zd			' Denominator
	
	Fraction = V		' Vrátíme pole o třech údajích
End Function

neutr · 10. 11. 2014 11:33:36

Tak se na to dívám a musím poznamenat, že to jde mnohem lépe. Jen pár postřehů :
----------------------------------------------------------------------------------
- Když iterujeme For i = 1 To N tak například N/i má ekvivalent v i/N. Postačí tedy iterovat do poloviny N s tím, že číslo > 1 je v normálním poměru N/i a číslo < 1 je v rozsahu i/N. Pokud se iteruje na polovny, je potřeba inkrementovat po dvojkách (1/2 ze 6 = 3 a stejně 1/2 ze 7 = 3,5 ale iterátor (For) si s tím poradí po svém - je to zase jen 3) - celočíselné hodnoty - takže lichá N přeskakovat. Pokud to uděláme mimo cykllus For tak máme zaděláno na problém. Při první liché polovině (například ze 7 = 3,5 a další lichá polovina už to posune na X,75.
----------------------------------------------------------------------------------
- Když budeme provádět přibližování (nejen pomocí PN) tak pomocí While-Wend - ale bez kontroly zda poměr vyhovuje. Stačí udělat INT(testované číslo)+1. Respektive pro malá čísla bychom mohli udělat
INT(testované číslo)+2. Pak
While (n-iVar)/iVar > INT(testované číslo)+X.
iVar +1 (+x)
Wend
Následně start s testem
For i = iVar To INT(N/2) Step + 2
a pak testovat podmínky.
....
....
Next i
----------------------------------------------------------------------------------
- Podmínky pro přiblížení s tolerancí - například 5% - zadáme jen y = 5.
Lim1 = ((testované číslo) / 100)*(100-y)
Lim2 = ((testované číslo) / 100)*(100+y)
Target = [ Lim1..Lim2 ] (ať už jako N/i, nebo jako desetinné číslo.
Testuje na poslední hodnotu "před Targetem" a také "za Targetem". Jakmile testovaný poměr skončí v "Targetu" můžeme skončit. K tomu nám postačují 2 proměnně T1 a T2 (T jako target).
Ovšem pokud nás bude například i někdy později zajímat přesný poměr, tak se vyplatí otestovat všechny hodnoty v oblasti "Target". Tento postup má následující filozofii :
T1, T2 = limita testovaného rozsahu - poměrně velká možná i desítky procent.
p(1), ..p(x) = testované hodnoty (poměry) jako prvky "n".
R = result jako co nejvíce přesné vyjádření poměru (limity stroje ap).
1. stav ...................... p(x)........T1.........R.......T2.......p(x+1)
2. n*2........................ T1..........p(x).......R.......p(x+1)...T2
3. T1=p(x), T2 = p(x+1), n*2.. T1..........p(x).......R.......p(x+1)...T2
4. T1=p(x), T2 = p(x+1), n*2.. T1..........p(x).......R.......p(x+1)...T2
A tak dál. To je přibližování se binární řadou. Dčíve jsem uvedl, že i inkremetace by měla být po dvou - a nyní ještě také to, že "n" má růst mocninou (násobkem) dvojky. Doufám, že je filozofie jasná:
----------------------------------------------------------------------------------
- Vždy mezi dvěma nějakými sousedními poměry N/i (respektive i/N) se vyskytne testovaný poměr. Tím, že postupujeme násobkem 2 (lépe binární) řadou, dělíme "Target" na polovinu (i když je to hyperbolická funkce, tak je to jedno). Znamená to, že pokud od samého začátku budeme testovat jen n*2 neuděláme chybu (ale také se nejspíš netrefíme do skutečného poměru). Ovšem existuje lepší cesta - zejména pro testovaný čísla < 1. Tou cestou je testování pomocí 2^n.
To má ještě jeden poměrně fatální důsledek. Nemusíme dělat žádné odhady a urychlení startovních hodnot. Najdeme - hledaný poměr mezi dvěma binárními výrazy a prvné před měl pořadí 22, nová řada (n*2) začne na 44 (na relativně stejném úseku je dvakrát více dílů). Pak také známe konec iterace - je to 46 - tedy For i = 44 to 46. - Testujeme jenom 2 úseky 44-45 a 45-46 (stačí jeden test - když ne první, pak druhý úsek). A tak dál.
Takže žádné složité věci na pochopení, ale něco také navíc. - Dokud se vyskytuje hledaní hodnota ve stejném pořadí poloviny "Targetu" - jde o sudé číslo. Jakmile se vyskytne v tom jiném dílu (polovině) - bude hrát v přesném poměru úlohu liché číslo. Pokud dojde někdy později k přechodu do původní poloviny - tak je zde další "sudý faktor" atd.
Z toho vlastně vyplývá jednoduché rozhodnutí - tolerance by měla být určena jako díl 2^n. Stejně tak startovní množina "N" - otestujeme INT(INT*(1/x) + (1/x)) - to je zaručeně menší, nežli skutečné, ale je blízko. Najdeme nejblíže vyšší 2^n. Chceme - li toleranci cca 6% postačí zřejmě startovat na "n" = 16 - najít target - určit konstantní start cyklu For. Použijeme vyšší "n" ze dvou možných . buď z tolerance, nebo z odhadu 2^n = nejblíže vyšší k "N"(z výrazu INT(...))
----------------------------------------------------------------------------------
Pro přesné dosažení poměru je nutné testovat i liché "N". Abychom se nezacyklovali, potřebujeme přesnější testování. Tam už musíme asi šetřit poměry mezi T1:R:T2. Odhadovat funkci a "N" řídit výpočtem.
----------------------------------------------------------------------------------
Takže testovat na limity "posouvaných" T1 a T2 z binárního vyjádření tolerance // startu 2^n nějblíže k "N". To nemůže trvat dlouho. Jiná je otázka nalezení skutečného poměru bez tolerance.
Moýná se do toho pustím, i když aby to bylo dobré, chce to hodně času zejména na testy. Moc se mi nelíbí While-Wend i když je rychlé. A tak je potřeba testovat i nějdivočejší čísla. Dokonce mám nepčíjemný dojem, že testovat jen jediný podíl na jednom šetřením čísle nemusí být správné. Představte si desetinné číslo jako součet, součin, podíl... ap. dvou či více jiných základních "zlomků". Viz například vyjádření Pí jako součtu.

Editoval neutr (10. 11. 2014 11:54:56)

OpenOffice.cz - web věnovaný kancelářským balíkům OpenOffice.org

#1 9. 11. 2014 10:30:52

Převod iracionálního čísla na zlomek - VYŘEŠENO

#2 9. 11. 2014 12:43:03

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

#3 9. 11. 2014 13:08:16

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

#4 9. 11. 2014 13:27:10

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

#5 9. 11. 2014 13:46:48

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

#6 9. 11. 2014 16:25:30

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

#7 10. 11. 2014 11:33:36

Re: Převod iracionálního čísla na zlomek - VYŘEŠENO

Zápatí