2012-12-14 7 views
7

Supponiamo di avere un numero reale. Voglio approssimarlo con qualcosa della forma a + sqrt (b) per i numeri interi aeb. Ma non conosco i valori di aeb. Ovviamente preferirei ottenere una buona approssimazione con piccoli valori di aeb. Lasciamo per ora indefinito cosa si intende per "buono" e "piccolo". Qualsiasi definizione sensata di questi termini andrà bene.convertire il numero reale in radicali

C'è un modo sano per trovarli? Qualcosa come il continued fraction algorithm per trovare approssimazioni frazionarie di decimali. Per ulteriori informazioni sul problema delle frazioni, vedere here.

MODIFICA: per chiarire, è un numero reale arbitrario. Tutto quello che ho sono un mucchio delle sue cifre. Quindi, a seconda dell'ottima approssimazione che vogliamo, a e b potrebbero o potrebbero non esistere. La forza bruta non è naturalmente un algoritmo particolarmente buono. Il meglio che posso pensare sarebbe di iniziare ad aggiungere numeri interi al mio reale, quadrando il risultato e vedendo se mi avvicino a un numero intero. Piuttosto forza bruta, e non un algoritmo particolarmente buono. Ma se non esiste nulla di meglio, sarebbe interessante da sapere.

MODIFICA: Ovviamente b deve essere zero o positivo. Ma un qualsiasi potrebbe essere intero.

+0

Due cose non sono chiare per me: 1- in quale forma hai inizialmente il tuo numero reale; 2- considerando che l'algoritmo naive a forza bruta fornisce la risposta migliore nel tempo O (r) assumendo operazioni aritmetiche a tempo costante, come si chiama "sano"? (a meno che tu non permetta al parametro a di essere negativo, nel qual caso la forza bruta non funziona) –

+0

Sì, sarebbe utile una migliore descrizione dei tuoi numeri reali. Alcuni real non possono essere rappresentati da un + sqrt (b), come Pi. – ckb

+0

A e b dovrebbero essere interi positivi? –

risposta

6

Nessuna necessità di frazioni continue; basta calcolare la radice quadrata di tutti i valori "piccoli" di b (fino a qualsiasi valore si ritenga ancora abbastanza "piccolo"), rimuovere tutto prima del punto decimale e ordinare/memorizzarli tutti (insieme allo b che lo ha generato).

Quindi quando è necessario approssimare un numero reale, trovare il radicale la cui porzione decimale è chiusa alla porzione decimale del numero reale. Questo ti dà b - scegliere il corretto a è quindi una semplice questione di sottrazione.

+0

Sono divertito come quando scrivevo la mia risposta ero principalmente interessato all'ottimizzazione per "buono", mentre si sta ottimizzando l'aspetto "piccolo". Non ho considerato questo approccio affatto. – amulware

0

Non so se esiste un qualche tipo di algoritmo standard per questo tipo di problema, ma mi incuriosisce, quindi ecco il mio tentativo di sviluppare un algoritmo che trovi l'approssimazione necessaria.

Chiama il numero reale in questione r. Quindi, per prima cosa suppongo che a possa essere negativo, in tal caso possiamo ridurre il problema e ora dobbiamo solo trovare un b in modo che la parte decimale di sqrt(b) sia una buona approssimazione della parte decimale di r. Scriviamo ora r come r = x.y con x essendo il numero intero e la parte decimale.

Now: 
b = r^2 
    = (x.y)^2 
    = (x + .y)^2 
    = x^2 + 2 * x * .y + .y^2 
    = 2 * x * .y + .y^2 (mod 1) 

Ora abbiamo solo per trovare un x tale che 0 = .y^2 + 2 * x * .y (mod 1) (circa).

Riempimento che x nelle formule sopra otteniamo b e possiamo quindi calcolare a come a = r - b. (Tutti questi calcoli devono essere accuratamente arrotondati.)

Ora, per il momento non sono sicuro se c'è un modo per trovare questo x senza forzarlo. Ma anche in questo caso, si può semplicemente usare un semplice ciclo per trovare uno x abbastanza buono.

Sto pensando a qualcosa di simile (pseudo codice semi):

max_diff_low = 0.01 // arbitrary accuracy 
max_diff_high = 1 - max_diff_low 
y = r % 1 
v = y^2 
addend = 2 * y 
x = 0 
while (v < max_diff_high && v > max_diff_low) 
    x++; 
    v = (v + addend) % 1 
c = (x + y)^2 
b = round(c) 
a = round(r - c) 

Ora, credo che questo algoritmo è abbastanza efficiente, mentre anche che consente di specificare la precisione desiderata della approssimazione. Una cosa che potrebbe essere fatta per trasformarla in un algoritmo O (1) è calcolare tutti gli x e inserirli in una tabella di ricerca. Se ci si preoccupa solo delle prime tre cifre decimali di r (ad esempio), la tabella di ricerca avrà solo 1000 valori, che sono solo 4kb di memoria (supponendo che siano utilizzati gli interi a 32 bit).

Spero che questo sia utile a tutti. Se qualcuno trova qualcosa di sbagliato nell'algoritmo, per favore fammelo sapere in un commento e lo aggiusterò.

MODIFICA: Riflettendo, ritiro la mia richiesta di efficienza. In effetti, non posso garantire che l'algoritmo sopra descritto finisca mai e, anche se dovesse accadere, potrebbe richiedere molto tempo per trovare un valore molto ampio di x che risolva adeguatamente l'equazione.

Si potrebbe forse tenere traccia dei migliori x trovati finora e rilassare i limiti di precisione nel tempo per assicurarsi che l'algoritmo si interrompa rapidamente, al possibile costo di precisione.

Questi problemi sono ovviamente inesistenti, se si calcola semplicemente una tabella di ricerca.

5

Questo è in realtà più un problema di matematica che un problema con il computer, ma per rispondere alla domanda penso che tu abbia ragione che puoi utilizzare le frazioni continue. Quello che fai è prima rappresentare il numero di destinazione come una frazione continua. Ad esempio, se si vuole approssimare pi (3.14159265) allora il CF è:

3: 7, 15, 1, 288, 1, 2, 1, 3, 1, 7, 4 ...

Il passo successivo è creare una tabella di CF per radici quadrate, quindi confrontare i valori nella tabella con la parte frazionaria del valore di destinazione (qui: 7, 15, 1, 288, 1, 2, 1, 3, 1 , 7, 4 ...). Ad esempio, supponiamo che il tuo tavolo abbia radici quadrate solo per 1-99. Quindi si troverà la corrispondenza più vicina sarebbe sqrt (51) che ha un CF di 7: 7,14 che si ripete. Il 7,14 è il più vicino a 7,15 di pi. Così la vostra risposta sarebbe:

sqrt (51) -4

Come la più vicina approssimazione dato un b < 100 che è fuori dal 0.00016. Se permetti b più grandi allora potresti ottenere una migliore approssimazione.

Il vantaggio dell'utilizzo di CF è che è più veloce di lavorare, diciamo, raddoppia o usa virgola mobile. Ad esempio, nel caso precedente devi solo confrontare due numeri interi (7 e 15), e puoi anche usare l'indicizzazione per rendere molto veloce la ricerca della voce più vicina nella tabella.

+0

Anche una risposta eccellente. Ha il vantaggio che può essere esteso ad altri problemi. Ma ho accettato la risposta che ritengo sia la migliore per la domanda, come affermato. –

1

Questo può essere fatto utilizzando mixed integer quadratic programming molto efficiente (anche se non ci sono garanzie run-time come MIQP è NP-completo.)

Definire:

d := the real number you wish to approximate 
b, a := two integers such that a + sqrt(b) is as "close" to d as possible 
r := (d - a)^2 - b, is the residual of the approximation 

L'obiettivo è quello di ridurre al minimo r. Imposta il tuo programma di quadratica come:

x := [ s b t ] 
D := | 1 0 0 | 
    | 0 0 0 | 
    | 0 0 0 | 
c := [0 -1 0]^T 
with the constraint that s - t = f (where f is the fractional part of d) 
and b,t are integers (s is not) 

Questa è una convessa (quindi in modo ottimale risolvibili) miscelato programma intero quadratica dal D è positiva semi-definita.

Dopo aver calcolato s,b,t, è possibile derivare la risposta utilizzando b=b, s=d-a e t può essere ignorata.

Il tuo problema potrebbe essere NP-completo, sarebbe interessante provare in tal caso.

1

Alcune delle risposte precedenti utilizzano metodi che sono di complessità temporale o spaziale O (n), dove n è il più grande "numero piccolo" che sarà accettato. Al contrario, il metodo seguente è O (sqrt (n)) nel tempo e O (1) nello spazio.

Supponiamo che numero reale positivo r = x + y, dove x=floor(r) e 0 ≤ y < 1. Vogliamo approssimare r per un numero del modulo a + √b. Se x+y ≈ a+√b quindi x+y-a ≈ √b, quindi √b ≈ h+y per alcuni offset intero h e b ≈ (h+y)^2. Per rendere b un intero, vogliamo minimizzare la parte frazionaria di (h+y)^2 su tutti i numeri idonei h. Ci sono al massimo √n valori ammissibili di h. Vedi il seguente codice Python e l'output di esempio.

import math, random 

def findb(y, rhi): 
    bestb = loerror = 1; 
    for r in range(2,rhi): 
     v = (r+y)**2 
     u = round(v) 
     err = abs(v-u) 
     if round(math.sqrt(u))**2 == u: continue 
     if err < loerror: 
      bestb, loerror = u, err 
    return bestb 

#random.seed(123456)  # set a seed if testing repetitively 
f = [math.pi-3] + sorted([random.random() for i in range(24)]) 
print (' frac  sqrt(b)  error  b') 
for frac in f:     
    b = findb(frac, 12) 
    r = math.sqrt(b) 
    t = math.modf(r)[0]   # Get fractional part of sqrt(b) 
    print ('{:9.5f} {:9.5f} {:11.7f} {:5.0f}'.format(frac, r, t-frac, b)) 

(Nota 1: Questo codice è in forma demo; i parametri per findb() sono y, la parte frazionaria di r, e rhi, la radice quadrata del numero più piccolo Si potrebbe desiderare di cambiare l'uso di. . parametri Nota 2: la linea
if round(math.sqrt(u))**2 == u: continue
di codice impedisce findb() di restituire valori perfetto-quadrati di b, ad eccezione del valore b = 1, perché nessun quadrato perfetto può migliorare l'accuratezza offerta da b = 1)

.

Segue l'output di esempio. Circa una dozzina di linee sono state elise nel mezzo. La prima riga di output mostra che questa procedura produce b=51 per rappresentare la parte frazionaria di pi, che è lo stesso valore riportato in alcune altre risposte.

frac  sqrt(b)  error  b 
    0.14159 7.14143 -0.0001642  51 
    0.11975 4.12311 0.0033593  17 
    0.12230 4.12311 0.0008085  17 
    0.22150 9.21954 -0.0019586  85 
    0.22681 11.22497 -0.0018377 126 
    0.25946 2.23607 -0.0233893  5 
    0.30024 5.29150 -0.0087362  28 
    0.36772 8.36660 -0.0011170  70 
    0.42452 8.42615 0.0016309  71 
    ... 
    0.93086 6.92820 -0.0026609  48 
    0.94677 8.94427 -0.0024960  80 
    0.96549 11.95826 -0.0072333 143 
    0.97693 11.95826 -0.0186723 143 

Con il seguente codice aggiunto alla fine del programma, appare anche l'output mostrato di seguito. Questo mostra approssimazioni più ravvicinate per la parte frazionaria di pi.

frac, rhi = math.pi-3, 16 
print (' frac  sqrt(b)   error   b  bMax') 
while rhi < 1000: 
    b = findb(frac, rhi) 
    r = math.sqrt(b) 
    t = math.modf(r)[0]   # Get fractional part of sqrt(b) 
    print ('{:11.7f} {:11.7f} {:13.9f} {:7.0f} {:7.0f}'.format(frac, r, t-frac, b,rhi**2)) 
    rhi = 3*rhi/2 

    frac  sqrt(b)   error   b  bMax 
    0.1415927 7.1414284 -0.000164225  51  256 
    0.1415927 7.1414284 -0.000164225  51  576 
    0.1415927 7.1414284 -0.000164225  51  1296 
    0.1415927 7.1414284 -0.000164225  51  2916 
    0.1415927 7.1414284 -0.000164225  51  6561 
    0.1415927 120.1415831 -0.000009511 14434 14641 
    0.1415927 120.1415831 -0.000009511 14434 32761 
    0.1415927 233.1415879 -0.000004772 54355 73441 
    0.1415927 346.1415895 -0.000003127 119814 164836 
    0.1415927 572.1415909 -0.000001786 327346 370881 
    0.1415927 911.1415916 -0.000001023 830179 833569