Quizz: Pire ailleurs?

Samuel Devulder · **Posté:** 10 Nov 2015, 08:30

Voici un nouveau quizz: que fait ce programme de 251 octets qui tourne aussi bien sur TO que MO?

Code:

****************************************
* debut  : $8000
* fin    : $80FA
* taille : 251
****************************************

        org     $8000

init    fcb     $34,$7F,$CC,$80,$02
        fcb     $1F,$8B,$CE,$1F,$40
        fcb     $A6,$C4,$63,$C4,$A1
        fcb     $C4,$26,$04,$86,$81
        fcb     $97,$BA,$A7,$C4,$4F
        fcb     $CE,$80,$FB,$ED,$C1
        fcb     $11,$83,$95,$0B,$26
        fcb     $F8,$8E,$03,$01,$4F
        fcb     $5F,$34,$56,$10,$8E
        fcb     $14,$11,$8E,$0A,$08
        fcb     $EC,$C3,$58,$49,$58
        fcb     $49,$E3,$C4,$58,$49
        fcb     $DD,$56,$34,$10,$A6
        fcb     $E4,$E6,$63,$3D,$D7
        fcb     $52,$EC,$61,$3D,$D7
        fcb     $54,$A6,$61,$E6,$63
        fcb     $3D,$8B,$00,$8B,$00
        fcb     $C3,$00,$00,$31,$3E
        fcb     $8D,$63,$ED,$62,$AF
        fcb     $C4,$35,$10,$30,$1F
        fcb     $26,$CC,$4F,$10,$8E
        fcb     $00,$0A,$8D,$52,$AF
        fcb     $C4,$C1,$09,$26,$04
        fcb     $0C,$F9,$20,$28,$C1
        fcb     $0A,$26,$06,$0C,$FA
        fcb     $CC,$30,$00,$8E,$86
        fcb     $39,$1E,$89,$10,$8E
        fcb     $80,$BA,$E7,$E2,$C6
        fcb     $30,$DB,$FA,$8D,$21
        fcb     $97,$FA,$E6,$E0,$0D
        fcb     $F9,$27,$06,$8D,$17
        fcb     $0A,$F9,$26,$FA,$35
        fcb     $56,$30,$1F,$26,$82
        fcb     $8E,$80,$DE,$E6,$80
        fcb     $27,$04,$8D,$0C,$20
        fcb     $F8,$35,$FF,$E7,$E2
        fcb     $C6,$07,$8D,$02,$E6
        fcb     $E0,$3F,$82,$7E,$E8
        fcb     $03,$34,$26,$8E,$00
        fcb     $11,$CC,$00,$00,$59
        fcb     $49,$A3,$62,$24,$02
        fcb     $E3,$62,$69,$61,$69
        fcb     $E4,$30,$1F,$26,$F0
        fcb     $ED,$62,$35,$16,$43
        fcb     $53,$39,$2E,$2E,$2E
        fcb     $20,$3A,$2D,$29,$20
        fcb     $28,$63,$29,$20,$50
        fcb     $55,$4C,$53,$20,$4E
        fcb     $6F,$76,$20,$32,$30
        fcb     $31,$35,$0A,$0D,$00
        fcb     $00

        end     init

Pour ceux qui n'ont pas envie de saisir le code, voici un SAP qui le lance automatiquement.

Prehisto · **Posté:** 10 Nov 2015, 17:50

La moitié d'un pipi :lol:

(c'est-à-dire environ 3.1416)

Samuel Devulder · **Posté:** 10 Nov 2015, 19:23

Oui

c'est pour ca que l'algo fait plic-ploc à chaque nouvelle goutte.. :voyons:

Le prog produit les premières décimales de Pi avec un algo de type "goutte à goutte" (et donc le plic-ploc ci-dessus

) qui sort les décimales unes à unes. Les anglais appellent ca un algorithme "spigot", en référence à un robinet (spigot) qui fuit.

Les décimales arrivent quasiment en temps constant ce qui fait que le plic-ploc est très régulier. On peut s'en servir de métronome. A la louche, le programme sort une décimale toutes les deux secondes. Mais parfois (lorsqu'il y a des 9 dans les décimales), le ploc est plus long et une grosse goute apparait avec tous les chiffres manquants.

D'ailleurs à ce sujet, si on a le courage d'attendre jusqu'au bout (35mins?), le programme se termine sur une grosse série de neufs. Comme si PI était un nombre fini: 3.14xxx49999999999... = 3.14xxx5 tout court! D'où le smiley final trahissant le certain désarroi du programme et qui justifie le titre du fil: le PI RAILLEUR. Ce PI possiblement fini se moque de nous! :bouffon:

Mais (question subsidiaire), est-ce que PI possiblement fini est une erreur dans le code, ou est-ce qu'on nous ment depuis des siècles avec les infinies décimales du nombre ? :bah:

Prehisto · **Posté:** 10 Nov 2015, 22:03

Samuel Devulder a écrit:

Mais (question subsidiaire), est-ce que PI possiblement fini est une erreur dans le code, ou est-ce qu'on nous ment depuis des siècles avec les infinies décimales du nombre ? :bah:

S'il faut répondre à la question, je dirais ni l'un ni l'autre. En fait, il existe tout une ribambelle d'algorithmes pour calculer la valeur de Pi, des plus simples aux plus complexes, des plus anciens aux plus récents, des plus approximatifs aux plus précis. Je ne sais pas quel algorithme tu as choisi pour ta démonstration mais tu n'as pas dû choisir le plus récent (plusieurs milliers de milliard de décimales), qui aurait impliqué plus de calculs et plus de mémoire et en tout cas aurait outrepassé la contrainte des 256 octets. Et malheureusement, tous ces algorithmes donnent un nombre PI fini. Donc, que tu trouves un nombre PI fini n'a théoriquement rien d'étonnant.

J'ai juste ?

Samuel Devulder · **Posté:** 10 Nov 2015, 22:31

Qu'il soit fini parce que l'algorithme tourne en temps fini, c'est certain.

Mais là c'est 4999... qui apparait avec des 9 qui se suivent. Or 0.xxxxx999999 n'est pas un réel mais le décimal (fini) 0.xxxxxx(x+1).

C'est curieux ce 4999999 qui apparait comme ca avant la 770e décimale.

Mon petit doigt me dit que si personne n'avait remarqué cette zone (indice! :voyons:

), c'est que l'algo produit les mauvais chiffres et donc qu'il est buggé. A propos d'algo, le voici:

Code:

***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO ou MO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $8000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    pshs   d,x,y,u,cc,dp

* init registre direct-page + valeur tab
       ldd    #(*<-8)*256+2
       tfr    a,dp

* Detection TO/MO
       ldu    #$1F40
       lda    ,u
       com    ,u
       cmpa   ,u
       bne    macmo
* pour TO, transforme "goto 2"
* en "cmpa #2", lequel ne modifie
* aucun registre
       lda    #$81
       sta    <PUTC
macmo  sta    ,u
       
***************************************
* initialisation du tableau avec la
* valeur 2
       clra
       ldu    #tab
set2   std    ,u++
       cmpu   #tabend
       bne    set2

***************************************
* Bouble externe: un tour par chiffre
       ldx    #N
boucJ  clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,u
       
* x = i, y = 2*i+1
       ldy    #LEN*2+1
       ldx    #LEN

***************************************
* boucle interne

* recuperation tab[i]
boucI  ldd    ,--u

* tab[i]*10
       lslb
       rola
       lslb
       rola
       addd   ,u
       lslb
       rola
       std    <mul3+1
       
* calcul i*q
       pshs   x
       lda    ,s
       ldb    3,s
       mul
       stb    <mul1+1
       ldd    1,s
       mul
       stb    <mul2+1
       lda    1,s
       ldb    3,s
       mul
mul1   adda   #0
mul2   adda   #0

* tab[i]*10 + i*q
mul3   addd   #0 

* y = 2*i-1      
       leay   -2,y
       bsr    div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
       stx    ,u

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -1,x
       bne    boucI

* CC=0 pour la div
       clra
       ldy    #10
       bsr    div
*  d = q/10 tab[1] = q % 10
       stx    ,u
      
***************************************
* if q==9,
*    ++nines
       cmpb   #9
       bne    q_ne_9
       inc    <nines
       bra    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 cmpb   #10
       bne    q_ne_10
       inc    <predig
       ldd    #256*'0
*      bra    print
* remplace par ldx #nnnn
       fcb    $8E   
      
***************************************
* q!=9 && q!=10
*    print predig puis nines fois '9
q_ne_10 
       lda    #'9

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       ldy    #PUTC
       stb    ,-s
       ldb    #'0
       addb   <predig
       bsr    PLOP
       sta    <predig
       
       ldb    ,s+
       tst    <nines
       beq    nextJ
print2 bsr    PLOP
       dec    <nines
       bne    print2
       
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,u
       leax   -1,x
       bne    boucJ    
       
***************************************
* trois points de suspension et smiley
       ldx    #railleur
       
***************************************
* affiche la chaine terminee par 0
* pointee par x    
puts   ldb    ,x+
       beq    exit
       bsr    PUTC
       bra    puts

***************************************
* sortie       
exit   puls   d,x,y,u,cc,dp,pc

***************************************
* fait plop et affiche une lettre
***************************************
PLOP   stb    ,-s
       ldb    #7
       bsr    PUTC
       ldb    ,s+
* passe au travers

***************************************
* affiche un char
***************************************
PUTC   goto   2
       jmp    $E803

***************************************
* Divise CC:D par Y. 
*
* En retour:
*     D=quotient 
*     X=reste
***************************************
div    pshs   d,y
       ldx    #17
       ldd    #0  * surtout preserver la retenue
div0   rolb
       rola
       subd   2,s
       bcc    div1
       addd   2,s
div1   rol    1,s
       rol    ,s
       leax   -1,x
       bne    div0
       std    2,s
       puls   d,x
       coma
       comb
       rts

railleur
       fcc    /... :-) (c) PULS Nov 2015/
       fcb    10,13    
* en sortie nines=0, donc on 
* economise le "fcb 0" de fin
* de chaine
*      fcb    0
nines  fcb    0
predig fcb    0

tab    rmb    LEN*2
tabend set    *

       end    ini

Le truc astucieux est la division de 17 bits par 16 bits. Les 17 bits sont nécessaires pour représenter la valeur "10*tab[i] + q*i" après le 700e tour de l'algo. Ce 17e bit est stocké dans la retenue, et une légère modification de la division présentée sur: viewtopic.php?p=1097#p1097 permet de la réaliser. En fait, j'ai pas mal galéré pour comprendre qu'il fallait 17 bits et implémenter cette division en conséquence. Au final la division sur 17bits est super simple. La modification par rapport aux 16bits est triviale, et elle n'est pas plus lente que pour du 16 bits.

Il existe une implémentation de cet algorithme pour ZX ici: http://www.pouet.net/prod.php?which=66568, mais elle n'a que 500 et quelques chiffres. Ici on en a 769. Le fait de ne travailler que jusquà 500 fait que 16 bits et pas 17 sont suffisants pour les calculs internes. En outre dans la version ZX la multiplication et la division y sont réalisés de façon naïve (boucle avec addition ou soustraction). Résultat leur algo est court, mais il est lent, très lent. Il est aussi limité à 500 et quelques décimales. Avec cette version thomson, les 17 bits sont suffisants pour atteindre 1050 décimales. Pour obtenir la 1051e il faut 18bits et là CC n'est plus suffisant: il faut réaliser une multiplication et une division sur 8+16=24bits. C'est assez facile à faire. Et une fois fait, avec l'arithmétique sur 24 bits on dépasse de loin plusieurs milliers de chiffres de PI. Le facteur limitant sera alors l'espace mémoire accessible et le temps de calcul.

En effet le temps de calcul est proportionnel au carré du nombre de décimales. Si on veut 2x plus de décimales, il faudra 4x plus de temps. Pour passer de 1000 chiffres (1h) à 100 000 il faudra alors 10 000x fois plus de temps. C'est très très long (plus d'un an).

Je pense que les algos basés sur les séries de Gauss sont plus rapides et ne nécessitent pas tout de suite une arithmétique 24bits. Mais ils ne sont pas réguliers: les chiffres arrivent de plus en plus vite vers la fin. Or j'aime bien la régularité de l'algorithme SPIGOT. Le beep qui est envoyé à chaque chiffre est comme un métronome: on "entend" l'algorithme tourner avec régularité.

A noter: sur MO5 l'affichage de chr$(7) ne produit pas de BEEP contrairement au TO7. Le programme est donc silencieux sur cette plateforme.

Samuel Devulder · **Posté:** 11 Nov 2015, 00:12

Voici un programme C qui permet d'étudier l'algorithme

Code:

/**
 * spigot.c
 *
 * (c) Samuel DEVULDER, Nov 2015
 */
#include <stdio.h>
#include <stdlib.h>

typedef unsigned long long num;

int bitmax(num v) {
   int i=0;
   while(v) {++i;v>>=1;}
   return i;
}

int main(int argc, char **argv) {
   int i,j,N,LEN;
   int nines=0, predigit=0;
   num *a, bitsA=0, bitsX=0, bitsQ=0;
   
   N = atoi(argv[1]);
   printf("N=%d\n",N);
   if(N==0) N=770;
   LEN = (N*10)/3+5;
   a = malloc(sizeof(*a)*(LEN+1)); if(a==NULL) exit(-1);
   
   for(i=LEN;--i;) a[i]=2;
   for(j=N; --j;) {
      num q = 0,x;
      for(i=LEN;i;--i) {
         bitsX |= (x = 10*a[i] + q*i);
         bitsA |= (a[i] = x % (2*i-1));
         bitsQ |= (q = x / (2*i-1));
      }
      a[1] = q%10; q=q/10;
      switch(q) {
         case 9: ++nines; break;
         case 10: printf(" %d", predigit+1); for(i=0;i<nines;++i) printf("0"); predigit=0;nines=0; break;
         default: printf(" %d", predigit);   for(i=0;i<nines;++i) printf("9"); predigit=q;nines=0; break;
      }
      fflush(stdout);
   }
   printf(" %d\n",predigit);
   free(a);
   
   printf("%d bits for %c\n", bitmax(bitsX), 'X');
   printf("%d bits for %c\n", bitmax(bitsA), 'A');
   printf("%d bits for %c\n", bitmax(bitsQ), 'Q');
   
   return 0;
}

Il affiche les chiffres de PI correspondant à l'entier passé en argument. Il affiche aussi le nombre de bits nécessaires pour X=tab[i]*10+q*i, pour A=les elements de tab[], et Q le quotient.

N=38, x=13bits, tab=8bits, q=7bits
N=533, x=16bits, tab=12bits, q=7bits (taille pour le ZX)
N=1050, x=17bits, tab=13bits, q=7bits
N=9830, x=21bits, tab=16bits, q=7bits

Plus généralement parlant le nombre de bits de x et des éléments du tableau augmente de 1 quand N est doublé:

Code:

N=512
16 bits for X
12 bits for A
7 bits for Q
user    0m0.077s

N=1024
17 bits for X
13 bits for A
7 bits for Q
user    0m0.202s

N=2048
18 bits for X
14 bits for A
7 bits for Q
user    0m0.702s

N=4096
19 bits for X
15 bits for A
7 bits for Q
user    0m2.671s

N=8192
20 bits for X
16 bits for A
7 bits for Q
user    0m10.640s

N=16384
21 bits for X
17 bits for A
7 bits for Q
user    0m40.671s

N=32768
22 bits for X
18 bits for A
7 bits for Q
user    2m40.030s

On voit aussi que le temps (user) est multiplié par 4 quand N est doublé.

Avec N=1024, il faut 17 bits pour x. Donc 17+7=24bits suffiraient pour calculer N=1024 * 2^7=131072 chiffres, mais cela prendrait (2^7)^2=16384 fois plus de temps que pour N=1024.

Prehisto · **Posté:** 11 Nov 2015, 10:44

Samuel Devulder a écrit:

A noter: sur MO5 l'affichage de chr$(7) ne produit pas de BEEP contrairement au TO7. Le programme est donc silencieux sur cette plateforme.

Même en changeant l'état du bit 3 de $2019 ?

Samuel Devulder · **Posté:** 11 Nov 2015, 11:49

ah oui c'était ca!

Concernant les 4999.. avec les 9 qui se poursuivent, c'est une séquence qui est parfaitement connue et on lui a donné le nom de point de Feynman°. Tant que je travaillais sur 16 bits je n'arrivais pas à obtenir cette séquence car les chiffres passés le 700e environ étaient faux. Le fait qu'on l'obtienne depuis que le programme stocke "tab[i]*10 + q*i" sur 17bits est la preuve que le programme est correct et qu'il marche bien.

Pour vraiment en être sur il faudrait envoyer le résultat une imprimante émulée et comparer le fichier-résultat avec les chiffres connus. Mais je suis confiant: en augmentant N jusqu'à 1050, on voit qu'après le point de Feynman les valeurs correspondent, alors que si l'algo échouait tous les chiffres suivants seraient faux à partir d'un moment. Avec N>1050 il faut 18bits et là l'algo produit les mauvais chiffres vers la fin.

Le programme est aussi 2x plus gros que la version ZX. Si on peut le réduire à 128 octets tout en restant aussi rapide (principalement division et multiplication), il pourrait être envoyé sur pouet pour faire un clin d'œil au ZX leur montrant que leur version n'est vraiment pas rapide

Samuel Devulder · **Posté:** 12 Nov 2015, 19:12

Héhé, j'ai posté un message sur pouet leur disant que c'était dommage de ne pas atteindre le point de Feynman, et du coup l'auteur vient de me dire que sisi ils on un autre algo, toujours en 128 octets qui fait les calculs sur 17 bits et qui dont atteint ce point. Il a l'air content du challenge.

Dans l'interim j'ai étudié leur algo qui n'est pas fondamentalement différent du miens. Il est plus compact parce que sur Zx il y a la possibilité d'additionner deux registres 16 bits directement, chose qu'on ne sait pas complètement faire sur 6809 (ajouter D à X est possible, mais ajouter X à Y ou X à X non). Et si on arrive à le faire, alors ca donne des instructions longues et lentes. Leur algo bénéficie aussi du fait qu'il y a deux jeux de registres parallèle sur ce CPU ce qui leur évite de stocker pleins de valeurs en mémoire. Toutes les instructions sont sur 1 seul octet.

Je réfléchis à la façon soit de rendre l'algo 6809 plus compact, soit plus rapide. Pour la compacité, c'est pas évident: il y a 4-5 octets gagnable vers la fin autour du cmpb #10, mais pas tellement plus. Or ca déborde de loin les 128 octets. Même en retirant le smiley et le support MO5 j'ai encore autour de 170 octets

Il n'y a pas à dire leur codage est super compact (car en plus ils affichent bien "3.14" au début et pas "0314").

Samuel Devulder · **Posté:** 12 Nov 2015, 23:24

Si je ne veux pas sacrifier la division sur 17bits, je descends à 160 octets et je ne vois pas où je peux grapiller les 32 octets en trop

Code:

(main)SMALLPI
***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $8000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    pshs   d,x,y,u,cc,dp

* init registre direct-page + valeur tab
       ldu    #predig
       pulu   d,dp,x,y

***************************************
* initialisation du tableau avec la
* valeur 2
set2   std    ,u++
       leax   -1,x
       bne    set2

***************************************
* Bouble externe: un tour par chiffre
       ldx    #N
boucJ  clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,y,u
       
* x = i, y = 2*i+1
       ldx    #LEN

***************************************
* boucle interne

* recuperation tab[i]
boucI  ldd    ,--u

* tab[i]*10
       lslb
       rola
       lslb
       rola
       addd   ,u
       lslb
       rola
* pile tab[i]*10,i,q
       pshs   d,x
* calcul i*q
       lda    2,s
       ldb    5,s
       mul
       stb    <mul1+1
       ldd    3,s
       mul
       stb    <mul2+1
       lda    3,s
       ldb    5,s
       mul
mul1   adda   #0
mul2   adda   #0

* tab[i]*10 + i*q
       addd   ,s++

* y = 2*i-1      
       leay   -2,y
       jsr    <div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
       stx    ,u

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -1,x
       bne    boucI

* CC=0 pour la div
       clra
       leay    10,x
       jsr    <div
*  d = q/10 tab[1] = q % 10
       stx    ,u
      
***************************************
* if q==9, 
*    ++nines
       inc    <nines
       addd   #$3930
       cmpb   #$39
       beq    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 blt    print
       inc    <predig
       ldd    #$3030
      
***************************************
* q!=9 && q!=10
*    print predig puis nines fois '9
q_lt_9 lda    #$39

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       stb    <print2+4
       ldb    <predig
       sta    <predig
print2 jsr    $E803
       ldb    #0
       dec    <nines
       bne    print2
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,y,u
       leax   -1,x
       bne    boucJ    

***************************************
* sortie       
exit   puls   d,x,y,u,cc,dp,pc

***************************************
* Divise CC:D par Y. 
*
* En retour:
*     D=quotient 
*     X=reste
***************************************
div    pshs   d,y
       ldx    #17
       ldd    #0  * surtout preserver la retenue
div0   rolb
       rola
       subd   2,s
       bcc    div1
       addd   2,s
div1   rol    1,s
       rol    ,s
       leax   -1,x
       bne    div0
       std    2,s
       puls   d,x
       coma
       comb
       rts

nines  fcb    0
predig fcb    0
* push order: CC D DP X Y U
regs   fcb    2
       fcb    *<-8
       fdb    LEN
       fdb    2*LEN+1
tab    rmb    LEN*2
tabend set    *

       end    ini

Prehisto · **Posté:** 12 Nov 2015, 23:48

A prime abord, je ne vois pas non plus. Si encore tu avais un ou deux octets à gagner, mais 32 ... :voyons:

Samuel Devulder · **Posté:** 13 Nov 2015, 01:16

Si je sacrifie le div par une implementation naïve ca n'est pas si lent que ca. C'est même plus rapide qu'avec la division classique.

Si je sacrifie le mul par une implementation lente, mais rikiki, j'arrive à un code bien plus lent mais plus petit.

Code:

(main)SMALLPI
***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $8000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    pshs   d,x,y,u,cc,dp

* init registre direct-page + valeur tab
       ldu    #predig
       pulu   d,dp,x,y

***************************************
* initialisation du tableau avec la
* valeur 2
*      ldd    #2
set2   std    ,u++
       leax   -1,x
       bne    set2

***************************************
* Bouble externe: un tour par chiffre
       ldx    #N
boucJ  clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,y,u
       
* x = i, y = 2*i+1
       ldx    #LEN

***************************************
* boucle interne

* recuperation tab[i]
boucI  ldd    ,--u

* tab[i]*10
       lslb
       rola
       lslb
       rola
       addd   ,u
       lslb
       rola
* pile tab[i]*10,i,q
       pshs   d,x
* calcul i*q
* lent mais court
       clra
       clrb
mul    addd   4,s
       leax   -1,x
       bne    mul
/
* rapide mais prend de la place
mul    lda    2,s
       ldb    5,s
       mul
       stb    <mul1+1
       ldd    3,s
       mul
       stb    <mul2+1
       lda    3,s
       ldb    5,s
       mul
mul1   adda   #0
mul2   adda   #0
/

* tab[i]*10 + i*q
       addd   ,s++

* y = 2*i-1      
       leay   -2,y
       jsr    <div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
       stx    ,u

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -1,x
       bne    boucI

* CC=0 pour la div
       clra
       leay    10,x
       jsr    <div
*  d = q/10 tab[1] = q % 10
       stx    ,u
      
***************************************
* if q==9, 
*    ++nines
       inc    <nines
       addd   #$3930
       cmpb   #$39
       beq    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 blt    print
       inc    <predig
       ldd    #$3030

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       stb    <print2+4
       ldb    <predig
       sta    <predig
print2 jsr    $E803
       ldb    #0
       dec    <nines
       bne    print2
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,y,u
       leax   -1,x
       bne    boucJ    

***************************************
* sortie       
exit   puls   d,x,y,u,cc,dp,pc

div    pshs   d,y
       ldx    #-1
       bcc    loop1
       bsr    loop2
loop1  bsr    loop2
       addd   2,s
       std    2,s
       stx    ,s
       puls   d,x,pc
loop2  leax   1,x
       subd   4,s
       bcc    loop2
       rts     
    
nines  fcb    0
predig
*      fcb    0
* push order: CC D DP X Y U
regs   fdb    2
       fcb    *<-8
       fdb    LEN
       fdb    2*LEN+1
tab    rmb    LEN*2
tabend set    *

       end    ini

Mais zut, c'est encore trop gros de 13 octets (141octets au total)

Bon là je ne vois plus quoi sacrifier. Peut-être qu'en laissant tomber les pshs et en utilisant un max de global...

Samuel Devulder · **Posté:** 13 Nov 2015, 11:40

Hum... je pourrais réduire encore un peu si je trouvais faire en sore que la carry résultant du calcul de tab[i]*10+q*i (tab[i]=--u, q=2,s i=x) par

Code:

       ldd    ,--u
* calcul i*q+tab[i]
mul    addd  2,s
       leax  -1,x
       bne   mul
* calcul i*q+tab[i]*10
       leax  9,x      
mul2   addd  ,u
       leax  -1,x
       bne   mul2      

était bonne. Le soucis c'est qu'on récupère la dernière carry de la dernière addition, or elle peut provenir de bien avant. Il faudrait avoir le moyen de la préserver une fois qu'elle apparait. J'ai bien pensé à encadrer les "addd" par un "push cc" et "orcc ,s+" qui auraient bien fait l'affaire... mais "orcc ,s+" n'existe pas.

[EDIT] trouvé: avec un ror je stock les carries du mul2. En sortie on a CC courant = la CC de la permière addition, et la zone mémoire des ROR = l'ensemble des 8 autres carries. Donc si la zone mémoire n'est pas nulle, il suffit de lever la carry en sortie. Tester si c'est différent de 0 et lever la carry peut se faire avec un simple neg. Youpie! :cool:

Résultat: 138 octets.. 10 de trop. Ca progresse

Code:

(main)SMALLPI
***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $8000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    pshs   d,x,y,u,cc,dp

* init registre direct-page + valeur tab
       ldu    #predig
       pulu   d,dp,x,y

***************************************
* initialisation du tableau avec la
* valeur 2
*      ldd    #2
set2   std    ,u++
       leax   -1,x
       bne    set2

***************************************
* Bouble externe: un tour par chiffre
       ldx    #N
boucJ  clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,y,u
       
* x = i, y = 2*i+1
       ldx    #LEN

***************************************
* boucle interne

* recuperation tab[i]
boucI  pshs   x
       ldd    ,--u
* calcul i*q+tab[i]
mul    addd  2,s
       leax  -1,x
       bne   mul
* calcul i*q+tab[i]*10
       leax  9,x      
mul2   addd  ,u
       ror   2,s
       leax  -1,x
       bne   mul2      
       neg  2,s

* y = 2*i-1      
       leay   -2,y
       jsr    <div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
       stx    ,u

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -1,x
       bne    boucI

* CC=0 pour la div
       clra
       leay    10,x
       jsr    <div
*  d = q/10 tab[1] = q % 10
       stx    ,u
      
***************************************
* if q==9, 
*    ++nines
       inc    <nines
       addd   #$3930
       cmpb   #$39
       beq    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 blt    print
       inc    <predig
       ldd    #$3030

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       stb    <print2+4
       ldb    <predig
       sta    <predig
print2 jsr    $E803
       ldb    #0
       dec    <nines
       bne    print2
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,y,u
       leax   -1,x
       bne    boucJ    

***************************************
* sortie       
exit   puls   d,x,y,u,cc,dp,pc

div    pshs   y
       ldx    #-1
       bcc    loop1
       bsr    loop2
loop1  bsr    loop2
       addd   ,s++
       exg    x,d
       rts
loop2  leax   1,x
       subd   2,s
       bcc    loop2
       rts     
    
nines  fcb    0
predig
*      fcb    0
* push order: CC D DP X Y U
regs   fdb    2
       fcb    *<-8
       fdb    LEN
       fdb    2*LEN+1
tab    rmb    LEN*2
tabend set    *

       end    ini

Samuel Devulder · **Posté:** 13 Nov 2015, 14:40

Ah, j'atteins 129 octets.. juste 1 de trop... on y est tout presque... suspens...

L'idée pour gagner est de préserver la valeur N dans X issu du LDU d'origine. Pour cela on déborde le tableau avec des $0002 jusqu'à ce que qu'on tape dans la ROM. Le débordement c'est pas gênant car l'algo travail en relatif sur [U-LEN*2, U[

J'ai aussi viré la sauvegarde du contexte en entrée et en sortie: le programme ne rends plus la main (de toute façon a on cassé le basic avec l'écriture en $9FFF-$DFFF).

Code:

(main)SMALLPI
***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $8000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    set    *
* init registre direct-page + valeur tab
       ldu    #regs
       pulu   d,dp,x,y

***************************************
* initialisation du tableau avec la
* valeur 2

* D,X: déjà chargés par le pulu
* U pointe du coup sur tab
set2   std    ,u++
       cmpd   -2,u
       beq    set2
* on déborde de loin la tab et on
* s'arrete quand on touche la rom.
* C'est pas très grave car le but
* est d'avoir LEN*2 octets de libres.
* Mais du coup X peut conserver sa
* valeur N, ce qui fait gagner 3 octets.

***************************************
* Bouble externe: un tour par chiffre
*       ldx    #N
boucJ   set    *
* clra inutile car la dernière carry
* stocké ne dépasse pas 10
*      clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,y,u
       
* x = 2*i, y = 2*i+1
       leax   -1,y

***************************************
* boucle interne

* recuperation tab[i]
boucI  pshs   x
       ldd    ,--u
* calcul i*q+tab[i]
mul    addd  2,s
       leax  -2,x
       bne   mul
* calcul i*q+tab[i]*9
       leax  9,x      
mul2   addd  ,u
       ror   2,s
       leax  -1,x
       bne   mul2      
       neg  2,s

* y = 2*i-1      
       leay   -2,y
       jsr    <div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
* stocké sur 0-2,u par le div

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -2,x
       bne    boucI

* CC=0 pour la div
       clra
       leay    10,x
       jsr    <div
*  d = q/10 tab[1] = q % 10
* stocké sur 0-2,u par le div
      
***************************************
* if q==9, 
*    ++nines
       inc    <nines
       addd   #$3930
       cmpb   #$39
       beq    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 blt    print
       inc    <predig
       ldd    #$3030

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       stb    <print2+4
       ldb    <predig
       sta    <predig
print2 jsr    $E803
       ldb    #0
       dec    <nines
       bne    print2
     
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,y,u
       leax   -1,x
       bne    boucJ    

***************************************
* sortie       
exit   bra    exit

* divise CC:D par Y
* d     : quotient
* 0-2,u : reste
div    pshs   y
       ldx    #-1
       bcc    loop1
       bsr    loop2
loop1  bsr    loop2
       addd   ,s++
       std    ,u
       tfr    x,d
       rts
loop2  leax   1,x
       subd   2,s
       bcc    loop2
       rts     
    
* push order: (CC) D DP X Y (U)
regs   fdb    2
       fcb    *<-8
       fdb    N
       fdb    2*LEN+1

nines  equ    regs
predig equ    regs+1
       
tab    rmb    LEN*2
tabend set    *

       end    ini

Quand je compare la version 6809 avec la toute dernière Z80, je vois qu'on est tombé sur quasiment la même division 17 bits simplifiée. Mais la leur de fait pas de "ADDD 2,s" à la fin pour contrebalancer le fait qu'on a fait une soustraction de trop quand la carry passe à 1. C'est bizarre.

Code:

; Division 17 bits Z80
L_A073  LD   BC,#FFFF
        CALL C,L_A079
L_A079  AND   A
L_A07A  INC   BC
        SBC   HL,DE
        JR   NC,L_A07A
        RET

Samuel Devulder · **Posté:** 13 Nov 2015, 16:17

Pour gagner l'octet restant il faudrait remplacer le "exit bra exit" par un truc en un seul opcode. Je pense au fameux HCF, mais je ne trouve plus son opcode. C'est lequel ??

[EDIT] Ah non plus besoin. En remplacant le cmpd -2,u par un simple cmpb -3,u qui fait la même chose (savoir si on a tapé dans la rom) j'obtiens

Code:

****************************************
* debut  : $C000
* fin    : $C07F
* taille : 128
****************************************

        org     $C000

init    fcb     $CE,$C0,$79,$37,$3E
        fcb     $ED,$C1,$E1,$5F,$27
        fcb     $FA,$5F,$34,$76,$30
        fcb     $3F,$34,$10,$EC,$C3
        fcb     $E3,$62,$30,$1E,$26
        fcb     $FA,$30,$09,$E3,$C4
        fcb     $66,$62,$30,$1F,$26
        fcb     $F8,$60,$62,$31,$3E
        fcb     $9D,$60,$ED,$62,$35
        fcb     $10,$30,$1E,$26,$DE
        fcb     $4F,$31,$0A,$9D,$60
        fcb     $0C,$79,$C3,$39,$30
        fcb     $C1,$39,$27,$18,$2D
        fcb     $05,$0C,$7A,$CC,$30
        fcb     $30,$1E,$89,$D7,$53
        fcb     $D6,$7A,$97,$7A,$BD
        fcb     $E8,$03,$C6,$00,$0A
        fcb     $79,$26,$F7,$35,$76
        fcb     $30,$1F,$26,$AD,$20
        fcb     $FE,$34,$20,$8E,$FF
        fcb     $FF,$24,$02,$8D,$09
        fcb     $8D,$07,$E3,$E1,$ED
        fcb     $C4,$1F,$10,$39,$30
        fcb     $01,$A3,$62,$24,$FA
        fcb     $39,$00,$02,$C0,$03
        fcb     $01,$14,$11

        end     init

YOUPIE YOUPIE YOUPIE :bien:

Voici le code final:

Code:

(main)SMALLPI
***************************************
* Calcul de PI avec la methode SPIGOT
* http://tinyurl.com/spigot-pdf
*
* Inspiré de:
* http://rosettacode.org/wiki/Pi#Lua
*
* C'est un algo "goute à goute"
* (spigot), qui délivre les décimales
* unes à unes.
*
* Fonctionne sur TO.
*
* (c) Samuel DEVULDER, Nov 2015
*
***************************************

       org    $C000

* nombre de chiffres
N      set    769

* taille du tableau a utiliser
* (5 valeurs en rab pour la precision)
LEN    set    (10*N/3)+5

       setdp  *<-8

***************************************
* Point d'entree       
ini    set    *
* init registre direct-page + valeur tab
       ldu    #regs
       pulu   d,dp,x,y

***************************************
* initialisation du tableau avec la
* valeur 2

* D,X: déjà chargés par le pulu
* U pointe du coup sur tab
set2   std    ,u++
       cmpb   -1,u
       beq    set2
* on déborde de loin la tab et on
* s'arrete quand on touche la rom.
* C'est pas très grave car le but
* est d'avoir LEN*2 octets de libres.
* Mais du coup X peut conserver sa
* valeur N, ce qui fait gagner 3 octets.

***************************************
* Bouble externe: un tour par chiffre
*       ldx    #N
boucJ  
* clra
       clrb
* 0-1,s represente le quotient courant
* (initialise a zero)
       pshs   d,x,y,u
       
* x = 2*i, y = 2*i+1
       leax   -1,y

***************************************
* boucle interne

* recuperation tab[i]
boucI  pshs   x
       ldd    ,--u
* calcul i*q+tab[i]
mul    addd  2,s
       leax  -2,x
       bne   mul
* calcul i*q+tab[i]*9
       leax  9,x      
mul2   addd  ,u
       ror   2,s
       leax  -1,x
       bne   mul2      
       neg  2,s

* y = 2*i-1      
       leay   -2,y
       jsr    <div
* q = (tab[i]*10 + i*q) / (2*i-1)
       std    2,s
* tab[i] = (tab[i]*10 + i*q) % (2*i-1)
* stocké sur 0-2,u par le div

* iteration sur l'indice precedent
* du tableau
       puls   x
       leax   -2,x
       bne    boucI

* CC=0 pour la div
       clra
       leay    10,x
       jsr    <div
*  d = q/10 tab[1] = q % 10
* stocké sur 0-2,u par le div
      
***************************************
* if q==9, 
*    ++nines
       inc    <nines
       addd   #$3930
       cmpb   #$39
       beq    nextJ

***************************************
* if q==10, 
*   print predig+1 puis nines fois '0
q_ne_9 blt    print
       inc    <predig
       ldd    #$3030

***************************************
* affiche predig, puis le remplace
* par b puis affiche nines fois la
* lettre contenue dans le reg a.
print  exg    a,b
       stb    <print2+4
       ldb    <predig
       sta    <predig
print2 jsr    $E803
       ldb    #0
       dec    <nines
       bne    print2
     
***************************************
* fin de boucle sur les chiffres
nextJ  puls   d,x,y,u
       leax   -1,x
       bne    boucJ    

***************************************
* sortie : pas de sortie      
exit   bra    exit

* divise CC:D par Y
* d     : quotient
* 0-2,u : reste
div    pshs   y
       ldx    #-1
       bcc    loop1
       bsr    loop2
loop1  bsr    loop2
       addd   ,s++
       std    ,u
       tfr    x,d
       rts
loop2  leax   1,x
       subd   2,s
       bcc    loop2
       rts     
    
* push order: (CC) D DP X Y (U)
regs   fdb    2
       fcb    *<-8
       fdb    N
       fdb    2*LEN+1

nines  equ    regs
predig equ    regs+1
       
tab    rmb    LEN*2
tabend set    *

       end    ini

Attention par rapport à la toute première version, il est extrêmement lent car il n'utilise pas l'instruction "mul", mais du coup il tient dans 128 octets!!

Logicielsmoto.com

Quizz: Pire ailleurs?

Qui est en ligne