Sida 1 av 1

bästa heltalsdivisionen?

Postat: 13 februari 2010, 03:27:27
av jesse
Har suttir och räknat bra länge nu, men orkar inte räkna ut den mest optimala lösningen.

Jag ska beräkna en kvot, och här är förutsättningarna:

resultatet är 16-bitars unsigned integer.
täljare och nämnare är båda 32-bitars (nämnaren är signed, men den kommer att vara positiv).

täljare A är ett tal i storleksordningen 2^26 till 2^31
nämnare B är ett tal i storleksordningen 2^20 till 2^25

Det inbördes förhållandet mellan talen är A / (B * 64) = 1 plus/minus några procent.

kvoten jag vill få fram ska representeras av ett heltal, där 32768 motsvarar "1" i ekvationen ovan.

Matematiskt uttryckt:
K = A / B * (32768/64).

Exempel: A =335E+6, B = 5.158E+6. Då är "kvoten" = A/(B*64) = 1.0148071
svaret ska då bli (32768 * 1.0148071) = 33253.


Problemet är: hur gör jag denna 32-bitars heltalsdivision för att få fram så bra noggrannhet som möjligt?

jag har funderat...

K = A / (B * 64) * 32768; // fungerar inte då A/(B*64) avrundas till 0 eller 1

K = A / B *512; // kanske så här helt enkelt ???
analys : A / B blir ca 64 - alltså 6 bitars noggrannhet. för dåligt!

K = A * 512 / B; // går inte , owerflow.

K = (A / 4) / ( B / 2048 + 0.5) ; // nu börjar det likna något...
men det blir inte 16 bitars noggrannhet:

int(B / 2048 +0.5) = 2519. ==> K = 33247 eller en avvikelse på 0.018% ( drygt 12 bitars noggrannhet). worst case, då A = 2^26 och B = ca 2^20 blir B/2048 = 512. ==> 0.2% eller 9 bitars noggrannhet.

hur hittar jag den optimala divisionen?

Re: bästa heltalsdivisionen?

Postat: 13 februari 2010, 03:51:18
av jesse
K = A / (B / 512 + 0.5);

+0.5 var ju intelligent.... funkar ju inte i heltalssammanhang... :wall:
Jag menar givetvis:

K = A / ( (B + 256) / 512);

Här får jag 11 bitars noggrannhet... tror jag kan nöja mig med det.

Om jag kör

K = A*2 / ( (B + 128) / 256);

så blir det 12 bitar... men då närmar jag mig gränsen för tillåtna max-värden.
om nu A är ett unsigned int och B är signed int, så kommer beräkningen att ske som signed, antar jag... eh. resultatet i alla fall. A * 2 kan i extrema fall få MSB som etta, då kan det tolkas som negativt, om inte jag först omvandlar B till unsigned:

K = A*2 / ( ( (uint32_t)B + 128) / 256);

kan antagligen inte bli bättre... då A(max) = 2^31-1

har jag tänkt helt fel?

Re: bästa heltalsdivisionen?

Postat: 13 februari 2010, 13:52:47
av sodjan
Kan du ge ett par exempel på hur "A" resp "B" kan se ut ?
Och också gärna vad A och B igentligen respresenterar. Var kommer de ifrån ?

> K = A / B * (32768/64).

Jag är inte helt med på varför du inte skriver K = A / B * 512 !?

> Det inbördes förhållandet mellan talen är A / (B * 64) = 1 plus/minus några procent

Och det du igentligen vill veta är hur mycket kvoten A/B*64 avviker från 1 ?

Sen det här med sign/unsigned. Om båda talen i alla fall alltid representerar
positiva tal så är det ju inte så intressant. Varför inte bara köra med båda som
postiva 32 bitars tal rakt av ?

Men, som ja sa i början, det riktigt intressanta här är ju inte hur man gör en division,
det finns det massor av lösningar på, utan vad det *igentligen* är du försöker göra.
D.v.s vad representerar "A" resp "B" i verkligheten, så att säga.

Re: bästa heltalsdivisionen?

Postat: 13 februari 2010, 17:00:33
av jesse
Just nu för mig så var divisionen det intressanta. Omständigheterna runt om varierar vid olika applikatoner, men den här typen av beräkning återkommer gång på gång, så jag hade lite funderingar kring vilka sätt som fungerar bra och inte. Det luriga är väl att man lätt förlorar en massa bitar vid heltalsdivision om man inte ser upp, det är därför jag väljer att använda 32-bitars heltal även om resultatet aldrig blir bättre än 14 bitars noggrannhet. Men det kanske finns enklare / bättre / annorlunda sätt att göra det på??? Eller så är det inte så mycket att välja på? Har man gjort det några gånger kanske det inte ens är en fråga utan en självklarhet hur man gör?

Jag vill få fram en kvot att multiplicera med för att kunna kalibrera förstärkningen på en ADC.

B är resultatet av en AD-omvandlare som kan ge resultatet 0 - 32767 inläst 1024 gånger där alla resultat summerats, och därför blir talet mellan 0 och 32767*1024. Teoretiska noggrannheten är max 14 bitar och jag vill inte förlora i noggrannhet pga divisionen.

B i det här fallet är en känd ström som appliceras vid kalibrering. 100% (eller 32768) motsvarar 250A. ADC:n ska kalibreras vid några olika strömmar, bl.a. 5% 10% och 20% av max.

A ska vara det förväntade värdet om "gain" hade varit idealt, dvs * 1.0000000, jag skiftade A så långt åt vänster som möjligt, därför är A 64 gånger större än det ideala värdet.

konstanten K = A/(B*64) används sedan vid mätning av en okänd ström: ström = mätvärde * K, men eftersom K inte är ett flyttal utan ett heltal, så väljer jag att multiplicera K med 32768 för att få upp till 15 bitars noggrannhet. Därför blir K = A/(B*64)*32768 eller A/B*512 och ADC resultatet beräknas således ström = mätvärde * K / 32768

exempelvis: vid strömmen 50A (20% av max) förväntas en ideal ADC ge resultatet 32768*0.20 = 6553.6, värdet läses in 1024 gånger: Då bakgrundsbruset skapar en oversamplingeffekt kommer slutsumman att komma i närheten av 6553.6 * 1024 = 6 710 886. Jag sätter därför A till 6 710 886 *64 = 4.2949673E+08. Detta är alltså inte ett inläst värde, utan ett teoretiskt värde som jag har som en konstant i programmet.

Antag att förstärkningen avviker med -1.2% då kommer jag att få B till ungefär 6 630 356.

K blir då A/(B*64) = 1.012 i decimal form. K representerat i binära formen 1:15 (en heltalsbit och 15 "binaler") motsvaras av 16-bitars heltalet 33166 (som man får från 1.012 * 2^15)

när väl konstanten K är beräknad, görs en multiplikation vid läsning av ADC för att få korrekt värde:

ström = K * inläst värde.

sodjan skrev: > K = A / B * (32768/64).

Jag är inte helt med på varför du inte skriver K = A / B * 512 !?
.
Jag försöker förtydliga de olika bestånsdelarna i uttrycket. /64 görs för att A är 64 ggr större än B. 32768 multiplicerar jag med för att K ska få sina "binaler" och inte bara bli 1.

Re: bästa heltalsdivisionen?

Postat: 13 februari 2010, 21:27:14
av ErikJ
Vi utgår alltså ifrån:

Kod: Markera allt

2^26 <= A < 2^31
2^20 <= B < 2^25
Jag har skrivit uträkningarna nedan stegvis för att få plats med kommentarer, med reservation för eventuella skrivfel och tankefel. :D Givietvis kan uträkningarna slås ihop till "one-liners" (utom de med while-satser, förståss).

Kod 1: Den här har du redan skrivit. (mellan 11 och 15 bitars precision)

Kod: Markera allt

    uint32_t A2 = A;       // < 2^31 ( >= 2^26 )
    uint32_t B2 = B / 512; // < 2^16 ( >= 2^11 )
    
    uint16_t K = (uint16_t)( A2 / B2 ); // 31-16 = 15 bitars precision i best case
                                        // 11 bitars i worst case (pga 11-bit nämnare)
Kod 2: Den här har du också skrivit. (mellan 12 och 15 bitars precision)

Kod: Markera allt

    uint32_t A2 = 2*A;       // < 2^32 ( >= 2^27 )
    uint32_t B2 = B / 256;   // < 2^17 ( >= 2^12 )
    
    uint16_t K = (uint16_t)( A2 / B2 ); // 32-17 = 15 bitars precision i best case
                                        // 12 bitars i worst case (pga 12-bit nämnare)
Kod 3: Om du tillåter att mellanuträkningar sker i större variabler så är det enkelt. Om inte uint48_t finns (mycket möjligt) kan den ersättas med uint64_t, vanligtvis kallad "unsigned long long", eller annan större än 48-bitars unsigned integer. Här får vi full precision (15 bitar) alltid.

Kod: Markera allt

    uint48_t A2 = (uint48_t)A * 32768; // < 2^46 ( >= 2^41 )
    uint32_t B2 = B * 64;              // < 2^31 ( >= 2^26 )
    
    uint16_t K = (uint16_t)( A2 / B2 ); // 46-31 = 15 bitars precision alltid
                                        // 26-15 = 11 bitar tillgodo i B2
Kod 4: Det är ju onödigt att multiplicera båda talen först, som i kod 3. Här multipliceras bara A.

Kod: Markera allt

    uint48_t A2 = (uint48_t)A * 512;   // < 2^40 ( >= 2^35 )
    uint32_t B2 = B;                   // < 2^25 ( >= 2^20 )
    
    uint16_t K = (uint16_t)( A2 / B2 ); // 40-25 = 15 bitars precision alltid
                                        // 20-15 = 5 bitar tillgodo i B2
Kod 5a: Kod 4 med Avrundnings-kompensering ( + ½ nämnare i täljaren ).

Kod: Markera allt

    uint48_t A2 = (uint48_t)A * 512;   // < 2^40 ( >= 2^35 )
    uint32_t B2 = B;                   // < 2^25 ( >= 2^20 )
    
    // Avrundningskompensering +B2/2, och vi har MÅNGA MSB lediga för att
    // undvika overflow
    uint16_t K = (uint16_t)( (A2+B2/2) / B2 ); // 40-25 = 15 bitars precision alltid
                                               // 20-15 = 5 bitar tillgodo i B2
Kod 5b: One-liner av kod 5a.

Kod: Markera allt

    uint16_t K = (uint16_t)( ((uint48_t)A*512 + (uint32_t)B/2) / (uint32_t)B );
Kod 6: Om du inte tillåter större än 32-bitars mellanlagring, kan du automatiskt justera förstärkningen av A & B så att du får maximal precision inom hela intervallen. OBS! Ingen avrundningskompensering!

Kod: Markera allt

    uint32_t A2 = A;
    uint32_t B2 = B;
    uint16_t B_delare = 512; // 32768/64

    while( A2 < 0x80000000 ) // Så länge vi har noll(or) först (MSB)
    {
        A2 *= 2;
        B_delare /= 2;
    }
    //  //  //  //  //  //  //  A är nu < 2^32 ( >= 2^31 )
    uint32_t B2 = B / B_delare;      // < 2^17 ( >= 2^16 )
    
    uint16_t K = (uint16_t)( A2 / B2 ); // 32-17 = 15 bitars precision alltid
                                        // 16-15 = 1 bit tillgodo i B2
Kod 7: Kod 6 med avrundningskompensering.
Eftersom A förstärktes maximalt i kod 6, så måste vi lämna en bit ledig för att undvika eventuell överflöde. Och vi hade en bit tillgodo i B2 i kod 6, som kan offras, så precisionen blir fortfarande 15 bitar.

Kod: Markera allt

    uint32_t A2 = A;
    uint32_t B2 = B;
    uint16_t B_delare = 512; // 32768/64

    while( A2 < 0x40000000 ) // Så länge vi har mer än 1 nolla först (MSB)
    {
        A2 *= 2;
        B_delare /= 2;
    }
    //  //  //  //  //  //  //  A är nu < 2^31 ( >= 2^30 )
    uint32_t B2 = B / B_delare;      // < 2^16 ( >= 2^15 )
    
    // Avrundningskompensering +B2/2, och vi har en hel MSB ledig för att
    // undvika overflow
    uint16_t K = (uint16_t)( (A2+B2/2) / B2 ); // 31-16 = 15 bitars precision alltid
                                               // 15-15 = 0 bitar tillgodo i B2
Sedan beror det på processorarkitekturen om man vill optimera bort några multiplikationer eller divisioner med skiftningar eller andra trix, för att öka exekveringshastigheten och/eller minska kodstorleken.

Har jag tänkt fel?
Har jag räknat fel ± 1 bit på precisionen?

Vilken processor/µC används?
Är det något som hindrar dig från att använda 48/64-bitars variabler?

(Oj, vad långt det blev :D)

EDIT:0x8000 --> 0x80000000, 0x4000 --> 0x40000000

Re: bästa heltalsdivisionen?

Postat: 14 februari 2010, 04:31:16
av jesse
loopen i kod6 belyser ju principen hur man optimerar en division så den underlättar ju tänkandet lite: Låt A och B vara så stora som möjligt är väl kontentan helt enkelt.
Är det något som hindrar dig från att använda 48/64-bitars variabler?
Nja, har inte testat i min kompilator, men vill alltid helst jobba så utrymmessnålt och snabbt som möjligt, och det gäller även här. Så alltihop är ju en fråga om kompromiss, eller rättare sagt; smart utnyttjande av resurserna utan att relevant data går förlorat. Därför ville jag knäcka nöten med en max 32 bitars division här. Loopen (kod6) har jag faktiskt använt mig av en gång (fast i assembler) i en liknande situation. Jag gjorde divisionen "själv" direkt i processorns register, så där ville jag inte heller ha mer än 32 bitar.

Så efter lite tänkande kommer jag nog att nöja mig med kod 2 i det här fallet. ADC:n ger ju inte fler bitar ändå, så allt annat vore overkill. Orsaken att jag startade tråden var att jag helt enkelt just då inte orkade reda ut hur jag skulle shifta de båda talen för att få precision, eller ens räkna rätt, men jag löste det ju på vägen så att säga.

Vanligtvis tar jag ett blankt A4 papper och sedan fyller med formler, sedan räknar jag några exempel för att se att jag inte räknat fel. Men den här gången valde jag att göra kalkylerna i en tråd i forumet. Kanske för att vädra lite hur jag tänker. Det verkar ju som ingen haft invändningar på själva tänket, så jag får väl godkänt då, antar jag.

Re: bästa heltalsdivisionen?

Postat: 14 februari 2010, 08:43:47
av Icecap
Om plats är ett problem kan man enkelt "dubbelboka" minne med "union", då kan man låta t.ex. 2 st 16 bit värden dela plats med en 32 bit, man får bara se till att inte använda dom samtidig.

Re: bästa heltalsdivisionen?

Postat: 14 februari 2010, 15:03:02
av jesse
Det är ju inte så att 64 bitar tar så värst mycket mer plats än 32.... det är koden:

Jag testade att lägga in en enda ändring på en rad i programmet:

före:
milliamp = curr.summa1 * C_MAX / 32768 / curr.antal;

efter:
milliamp = (uint64_t)curr.summa1 * C_MAX / 32768 / curr.antal;

Före : Program: 10458 bytes
Efter: Program: 14980 bytes

så de 32 extra bitarna krävde 4522 bytes extra flashminne, och ökade programmets storlek med 43% :shock:
Och då är optimeringen satt till -Os

Jag undviker alltså 64-bitars heltal, om det är möjligt.

Re: bästa heltalsdivisionen?

Postat: 14 februari 2010, 21:39:00
av dangraf
Jag är bara lite nyfiken på vad du egentligen försöker optimera?

rubriken säger att du vill hitta ett så effektivt sätt att göa en heltalsdivition vilket gäller själva kalibreringsdelen. Men när jag läser inläggen så verkar det som att det är resultatet som du vill ha så exakt som möjligt när man ska använda sin konstant.

Tekniken du tänkte använda är att dela upp mätområdet i mellan 0-100% i 10 delar. Beror detta på att insignalen inte är helt linjär? På vilket sätt är den ej helt linjär?

Sen funderar jag även lite på insignalen, du samplar 1024 värden vid kalibreringen för att få så exakt resultat som möjligt. Hur är det i normalfallet då du omvandlar din insignal till ditt linjäriserade värde?

Det jag försöker fiska efter är att om du måste ta såpass många värden som 1024 för att få ett stabilt värde under kalibreringen så kanske det skulle löna sig med någon form av filtrering när du väl omvandlar din signal för att få så bra noggrannhet som möjligt.

Bara några funderingar..

Re: bästa heltalsdivisionen?

Postat: 15 februari 2010, 02:20:45
av 4kTRB
Newton-Raphson divsion är en metod som är snabb.
Kanske du kan nyttja den metoden?
Tror du kan iterera tills du fått önskad noggrannhet.

Re: bästa heltalsdivisionen?

Postat: 15 februari 2010, 12:49:54
av jesse
Nja, det är en gammal vana jag har att ta kalibrerings-samples under en tidsperiod, för att undvika t.ex påverkan av 50 Hz brum eller annat. I det här fallet är det en halleffektsensor som "svävar" lite fram och tillbaks, så jag vill få fram ett medelvärde.

med "Den bästa divisionen" menade jag både kodeffektiv (utrymme i Flash) och precis (noggrannhet). Tiden det tar att beräkna konstanterna spelar ingen roll.

Jag delar upp kalibreringen i fyra olika nivåer (inte 10). Det beror på att jag läser av den analoga signalen med olika förstärkning beroende på amplituden. Därför vill jag kalibrera för varje förstärkarsteg för sig.

4ktrb: jag ska kolla upp det!

Re: bästa heltalsdivisionen?

Postat: 15 februari 2010, 13:30:35
av 4kTRB
En algoritm Goldschmidt finns beskriven på wiki förutom NR.
http://en.wikipedia.org/wiki/Division_(digital)

Re: bästa heltalsdivisionen?

Postat: 15 februari 2010, 16:22:45
av jesse
spännande länk! Jättekul. Jag kommer nog att testa några av metoderna när jag får tid.

Re: bästa heltalsdivisionen?

Postat: 15 februari 2010, 16:27:49
av 4kTRB
Det är ju avancerade saker som kräver vana med matematik och sånt.
Om du har kunskaperna så kan det helt klart vara värt ett försök.

Re: bästa heltalsdivisionen?

Postat: 27 februari 2010, 21:13:13
av Swech
En annan infallsvinkel är följande:
Du vill mäta ett värde som skall vara max 65535.
Differensen är inte så stor. max 1-2 % tal max +- 655

Ta då och addera ihop den uppmätta differensen (förväntat - uppmätt) istället för varje mätning,
Summan får plats med 20 bitar
Division med 1024 är ju enkelt därefter..

Swech