Un jour de 1971, M. Hanson a déclaré : « Aujourd'hui, je vais vous apprendre à obtenir la racine carrée d'un nombre. Un moyen plus rapide que le division longue façon dont je suis censé enseigner. Je vais vous apprendre avec un mot. Et vous n'oublierez jamais. Puis il écrit au tableau
Deviner
et le souligne. Lorsque les rires se sont calmés, il a expliqué comment diviser la supposition en nombre d'origine, moyenne, répétition. Il a probablement mentionné Newton Raphson .
La méthode de Newton pour trouver les racines ou les zéros d'une équation utilise la dérivée de l'équation pour accélérer le processus. La racine carrée de 10 ressemble à trouver une valeur de X qui fait X^2 == 10 ou X^2 – 10 == 0. La dérivée de X ^ 2 - 10 est 2X, et le cœur de la boucle pour trouver des valeurs successives de X est
Xnouveau = X - ( (X^2 - 10) / (2X) )
qui, si vous jouez avec, est exactement le deviner formule, et converge comme ceci, pour une estimation initiale de 2 :
2 - estimation initiale
3.5 – appliquer x = (x + 10/x)/2
3.17857142857143
3.16231942215088
3.16227766044414
3.16227766016838
3.16227766016838 – six itérations pour une réponse stable !
Voici un programme JSL à démontrer.
f = Expr( x ^ 2 - V );
step = Eval Expr( x - Expr( Name Expr( f ) ) / Expr( Derivative( Name Expr( f ), x ) ) );
V = 1e2; // find the sqrt of V
x = -8; // initial guess.
xnext = .; // set up loop
n = 0; // counter
While( 1,
xnext = step;
If( Abs( x - xnext ) < Abs( 1e-14 * x ),
Break()
);
x = xnext;
Write( "\!n", n, " ", x );
n = n + 1;
);
Write( "\!ndone, a root of " || Char( Name Expr( f ) ) || " for V = " || Char( v ) || " is ", x );
La première ligne définit une expression, f, qui aidera à trouver la racine carrée de V. La deuxième ligne, step, est l'itération de Newton-Raphson pour un cycle ; step apparaît dans la boucle while où il évalue la valeur X suivante. La fonction Dérivée fonctionne très bien sur des polynômes simples, donnant 2X dans cet exemple. Le NameExpr() obtient la formule de f sans l'évaluer et le Expr() fait exécuter le NameExpr() (pour obtenir la formule) quand EvalExpr() trouve le Expr(). Le pas est défini sur x - f(x)/f'(x).
À l'intérieur de la boucle while(1), les pas sont effectués, imprimés et comptés. Le test vérifie la différence entre les valeurs X actuelles et suivantes ; si la valeur est très petite à la taille actuelle de X, elle ne se rapprochera probablement pas et peut osciller entre deux nombres très proches, elle sort donc de la boucle while.
0 -10.25
1 -10.0030487804878
2 -10.0000004646115
3 -10
fait, une racine de x ^ 2 - V pour V = 100 est -10
N'oubliez pas qu'il existe à la fois une racine carrée positive et une racine carrée négative ! Les valeurs de départ déterminent à quelle racine mène la convergence. À partir de -8 dépassé jusqu'à -10,25, puis s'est installé à -10.
(Pourquoi divise-t-il par la dérivée, f' ? Intuitivement, lorsque la pente de la ligne f est proche de zéro, le point suivant devra être plus éloigné pour modifier le résultat Y de f, et la division par la pente le fait.)
Le JSL ci-dessus ne gère pas bien les requêtes impossibles. Si on lui demande de trouver sqrt (-100), le JSL boucle, à la recherche d'un nombre complexe qu'il n'est pas équipé pour trouver.
Nombres complexes
Dans l'index de script JSL, la recherche de complexe produit deux exemples différents à moitié cuits. L'un est basé sur les classes JSL et l'autre sur les espaces de noms et les matrices. Je pense que la version matricielle est plus proche d'une solution robuste, et probablement plus rapide. La version class est alléchante, mais j'ai eu du mal à l'utiliser (vitesse, problèmes de mémoire). Étant donné que JMP n'a pas de surcharges d'utilisateurs pour les opérateurs, ni l'un ni l'autre ne sera facile à utiliser.
Cette JSL est très similaire à la JSL précédente. Il y a un gros morceau de code pour gérer l'arithmétique des nombres complexes, la même formule, plus ou moins, et maintenant il gère aussi sqrt(-100).
ns = New Namespace("complex", // adapted from scripting index example
{
make = Function( {a, b}, a || b ),
add = Function( {x, y, t=[0 0], u=[0 0], v=[0 0], w=[0 0]}, // an example has 5 terms
if( Is Number( x ),x = x || 0);
if( Is Number( y ),y = y || 0);
x + y + t + u + v + w;
),
sub = Function( {x, y, t=[0 0], u=[0 0], v=[0 0], w=[0 0]}, // an example has 3 terms. this sub,add trick is good enough for a demo.
if( Is Number( x ),x = x || 0);
if( Is Number( y ),y = y || 0);
x - y - t - u - v - w;
),
mul = Function( {x, y},
if( Is Number( x ),x = x || 0);
if( Is Number( y ),y = y || 0);
local:z = J( 1, 2 );
local:z[1] = x[1] * y[1] - x[2] * y[2];
local:z[2] = x[1] * y[2] + x[2] * y[1];
local:z;
),
div = Function( {x, y},
if( Is Number( x ),x = x || 0);
if( Is Number( y ),y = y || 0);
local:z = J( 1, 2 );
local:d = (y[1] ^ 2 + y[2] ^ 2);
local:z[1] = (x[1] * y[1] + x[2] * y[2]) / local:d;
local:z[2] = (x[2] * y[1] - x[1] * y[2]) / local:d;
local:z;
),
// https://math.stackexchange.com/questions/1397437/formula-for-raising-a-complex-number-to-a-power
// If a,b are real, then a+ib=r(cosθ+i*sinθ) where r^2=a^2+b^2 and tanθ=b/a, and (a + i*b)^N = r^N(cos(Nθ) + i*sin(Nθ))
pow = function( {c, n},
if( Is Number( c ),c = c || 0);
local:a = c[1];
local:b = c[2];
local:theta = ATan( local:b, local:a );
local:r = Sqrt( local:a ^ 2 + local:b ^ 2 );
local:r ^ n * (Cos( n * local:theta ) || Sin( n * local:theta ));
),
abs = function( {x},
Sqrt( x[1] ^ 2 + x[2] ^ 2 )
),
toChar = Function( {x,w=20,d=14}, // rounding (to 14 places) is not a good idea if the numbers are small ( 1e-9 will lose lots of precision)
Char( round(x[1],d),w ) || "+" || Char( round(x[2],d),w ) || "i"
),
convertExpression = function( {e},
substitute( Name Expr( e ),
Expr( Add() ), Expr( complex:add() ),
Expr( Subtract() ), Expr( complex:sub() ),
Expr( Multiply() ), Expr( complex:mul() ),
Expr( Divide() ), Expr( complex:div() ),
Expr( Power() ), Expr( complex:pow() )
)
)
}
);Ci-dessus ajouté une fonction pow() pour x ^ 2 ou x ^ 3 etc où x est un nombre complexe, une fonction abs() qui devrait peut-être s'appeler magnitude() mais correspond à la version non complexe comme abs, un peu étrange - recherche de paramètres supplémentaires pour add() et sub(), et une fonction convertExpression(). Un grand bravo à stackexchange/Michael Hardy pour une fonction pow() facile à saisir !
La fonction convertExpression transforme une expression JSL avec les opérateurs +, -, *, / et ^ en une expression équivalente qui appelle les fonctions d'espace de noms complexes. Les opérateurs binaires JSL sont assez agressifs pour transformer ((a+b)+(c+d)) en add(a,b,c,d) qui n'est pas une fonction binaire... il a quatre arguments. C'est pourquoi il existe des arguments d'option supplémentaires sur les fonctions complex:add et :sub.
Le morceau suivant de JSL déplace la boucle dans une fonction, anticipant le développement futur :
findzero = function( {x},
{xnext, n},
xnext = .; // set up loop
n = 0; // counter
While( n < 200,
xnext = step;
If( complex:Abs( complex:sub( x, xnext ) ) < 1e-14 * complex:abs( x ),
Break()
);
x = xnext;
n = n + 1;
);
If( n == 200, x = [. .] );
Eval List( {x, n} );
);findzero prend un point de départ x et renvoie une liste avec la solution itérée x et le nombre d'itérations n.
Et le code pour appeler la fonction… la formule et la dérivée sont les mêmes pour les nombres complexes.
f = Expr( x ^ 2 - V );
step = Eval Expr( x - Expr( Name Expr( f ) ) / Expr( Derivative( Name Expr( f ), x ) ) );
step = complex:convertExpression( Name Expr( step ) );
V = complex:make( -100, 0 ); // find the sqrt of V==-100
x = complex:make( 0, 1 ); // initial guess. Does not like 0i for a start value.
{x, n} = findzero( x );
Write( "\!ndone, a root of " || Char( Name Expr( f ) ) || " for V = " || complex:toChar( v, 10, 5 ) || " is ", complex:toChar( x, 10, 5 ) );
Eval( complex:convertExpression( Name Expr( f ) ) );La 3ème ligne appelle la fonction convertExpression pour changer les opérateurs +, -, *, /, ^ des versions de nombre réel aux versions de nombre complexe. V est défini sur une valeur complexe et le point de départ initial est défini sur un point en dehors de la droite des nombres réels. Une fois terminé, le journal affiche
fait, une racine de x ^ 2 - V pour V = -100+0i est 0+10i
[0 0]
V=-100 + 0i est le nombre réel -100 (sans composante imaginaire.) 0 + 10i est la réponse imaginaire pure sans composante réelle ; 10i * 10i est 100*(i^2), ou -100, comme prévu. Les racines complexes apparaissent en paires conjuguées ; La méthode de Newton en a trouvé une basée sur le point de départ. L'autre est à -10i.
Modifions le code pour garder une trace du nombre d'itérations n pour un groupe de points dans le plan complexe et le graphique qui, peut-être comme un tracé de contour, et utilisons cubique au lieu de quadratique :
f = Expr( x ^ 3 - V );
step = Eval Expr( x - Expr( Name Expr( f ) ) / Expr( Derivative( Name Expr( f ), x ) ) );
step = complex:convertExpression( Name Expr( step ) );
V = complex:make( -2, 0 ); // find the sqrt of V==-2
dt = New Table( "complex",New Column( "xre" ),New Column( "yim" ),New Column( "iterations" ),New Column( "basinr" ),New Column( "basini" ));
stepsize=1/32;
for(xreal=-2,xreal<=2,xreal+=stepsize,
for(yimag=-2,yimag<=2,yimag+=stepsize,
x = complex:make( xreal, yimag ); // initial guess. Does not like 0i for a start value.
{x, n} = findzero( x );
dt<<ADDROWS></ADDROWS>Tous les points du plan complexe qui itèrent jusqu'à la même racine forment un bassin. Ci-dessus, le générateur de graphes combine les valeurs réelles et imaginaires en une seule valeur pour identifier une racine.
L'image montre les trois bassins en rouge, cyan et bleu foncé. Tous les points de départ rouges convergent vers la même valeur. Les tailles de point plus grandes mettent plus de temps à converger. Les points gris n'ont pas convergé après 200 itérations :
Three basins. Larger dots where the basins touch take longer to converge.Trois bassins. Les points plus grands où les bassins se touchent mettent plus de temps à converger.
Doubler la taille de l'image de 16 000 à 65 000 points prend environ 10 minutes avec l'implémentation des nombres complexes JSL. Je veux créer 10 000 bitmaps de points 4Kx2K, j'ai donc besoin d'une autre approche.
Amas Beowulf
UN Amas Beowulf est une pile de vieux ordinateurs, pour la plupart sans écran ni clavier, mais connectés sur un réseau local à un ordinateur de contrôle. Dans mon cas, l'ordinateur de contrôle utilise JSL pour dire au cluster quoi faire et pour rassembler les résultats. Le cluster basé sur Linux exécute des programmes Python créés par JMP. L'ordinateur Windows exécutant JMP utilise PLINK de RunProgram pour ouvrir une ou plusieurs sessions SSH sur chaque ordinateur et fournir le programme et les paramètres nécessaires. Le programme Python s'exécute, puis stocke les résultats sur un disque NFS, puis se termine. Lorsque JMP voit le RunProgram se terminer, il rassemble les résultats et démarre le RunProgram suivant. 10 000 RunPrograms plus tard, il y a suffisamment d'images pour faire une vidéo de cinq minutes.
Ethernet hub in center, below power-strip. Not enough mice and monitors and keyboards, SSH is best.Concentrateur Ethernet au centre, sous la multiprise. Pas assez de souris, de moniteurs et de claviers, SSH est le meilleur.
Python a un type de données de nombre complexe natif, donc une partie du JSL ci-dessus disparaît. Le JSL ressemble maintenant à ceci ; vous verrez du Python intégré sous forme de chaîne dans un instant.
F = Expr( // f and e parameters are driven by the frame animation time
f6 * xOld ^ e6 + f5 * xOld ^ e5 + f4 * xOld ^ e4 + f3 * xOld ^ e3 + f2 * xOld ^ e2 + f1 * xOld ^ e1 - 1
); // <<<<<<<<<<< supply equation here <<<<<<<<<<<<<<
scale = 2 / 1;
xresolution = 1920 * scale;
yresolution = 1080 * scale;
// keep aspect ratio by chopping shorter dimension...
xcenter = -1;
ycenter = 0;
delta = 3;
ytop = ycenter - delta * yresolution / xresolution;
xleft = xcenter - delta;
xright = xcenter + delta;
ybottom = ycenter + delta * yresolution / xresolution;
Le bit ci-dessus comprend l'équation et certains paramètres de taille vidéo. Le bit ci-dessous décrit le réseau. Hostkey, server, cpus, running, gigs sont des index dans la liste condat (données de connexion) des informations sur le cluster.
// named indexes for condat
hostkey = 1; // from a failed connection, this makes the connection work
server = 2; // the machine names
cpus = 3; // from the lscpu command, below
running = 4; // dynamic: how many are running on this cpu, right now
gigs = 5; // from the free command, below
// the avahi-daemon may need to be stop/start/status sometimes,
// unfortunately that is hard to do except at server keyboard if you can't connect...
Le problème avec le démon existe toujours, mais je sais maintenant qu'il annonce les machines avec un suffixe -2 (etc) et je peux le résoudre sans connecter d'écran.
condat = {
{"4e:c3:8a:3b:e7:d4:bd:0d:85:62:36:40:17:15:b4:a3", "dell790a.local", 4, 0, 22},
{"00:c2:37:09:dc:d7:cf:fb:33:e3:14:e3:d6:fb:39:be", "dell790b.local", 4, 0, 14},
{"54:c8:1a:6a:50:e3:51:02:ed:1b:f8:fa:bd:f4:e4:17", "dell7010a.local", 8, 0, 30},
{"66:0d:ca:3c:1f:11:eb:45:d4:92:d0:0b:ae:01:b8:6b", "dell7010b.local", 8, 0, 30},
{"94:84:2f:a5:56:2d:db:16:b5:80:a9:0d:87:07:ee:53", "dell9010.local", 8, 0, 30},
{"b2:b7:b3:c3:83:63:11:e7:56:e7:c9:c2:e6:fc:f0:c9", "thinkb.local", 8, 0, 30},
{"2b:ca:fc:45:f0:64:cd:ac:70:ce:93:73:e8:84:2d:75", "thinkc.local", 8, 0, 30},
{"62:b7:67:23:8d:86:fe:80:1a:62:ec:d8:41:1b:b4:1b", "thinkd.local", 8, 0, 30},
{"8e:54:ee:b0:25:67:b8:5c:4d:ad:5e:e9:71:87:fd:71", "netflix.local", 2, 0, 30} // host for vm
};Le JSL ci-dessous parcourt la liste de condat et convertit les noms tels que dell9010 en leur adresse IP en pointillés (ou équivalent IPV6). Il empêche les échecs de recherche ultérieurs ; le cluster a fonctionné pendant environ trois jours et je ne sais toujours pas pourquoi le bureau de ma femme fournit le service de recherche DNS… Avant de me souvenir de GetAddrInfo, j'ai utilisé le code ping commenté, qui est un peu plus maladroit. Ma maison est beaucoup plus heureuse avec IPV6 maintenant qu'elle ne l'était il y a quelques années. Parfois, l'adresse IPV6 appartient entre crochets.
// on my local network, the lookups for xxx.local fail randomly,
// and more often when busy. pre-lookup the dotted ip now.
//
// better way? Get Addr Info( "thinkd.local" )[3][4]
//
for( itry = 1, itry < 10, itry += 1,
nprocs = 0;
MAXCPU = 0;
good = 0;
for( iCondat = 1, iCondat <= N Items( condat ), iCondat += 1,
nprocs += condat[iCondat][cpus]; // nprocs will divide up the work (stepsize)
condat[iCondat][running] = 0;
MAXCPU = Max( MAXCPU, condat[iCondat][cpus] ); // MAXCPU loops needed later
if( Starts With( condat[iCondat][server], "[" ),
good += 1, // else
// txt = Run Program(
// Executable( "PING.EXE"/*path probably not needed*/ ),
// Options( {"-n 1", condat[iCondat][server]} ),
// ReadFunction( "text" )
// );
// ip = Regex( txt, "\[\[([^]]+)\]]\", "\0" ); // typ: [2603:6081:1a02:1200:8973:e7be:605c:f753]
try(
ip = "[" || Get Addr Info( condat[iCondat][server] )[3][4] || "]";
if( !Is Missing( ip ),
Write( "\!n", condat[iCondat][server] || " --> " || ip );
condat[iCondat][server] = ip;
good += 1;
, //
Write( "\!n", condat[iCondat][server] || " lookup failed ", Left( txt, 70 ) )
); //
, // catch
show(exception_msg,iCondat,condat[iCondat][server],Get Addr Info( condat[iCondat][server] ));
);
);
);
if( good == N Items( condat ),
Break(),
Wait( .1 )
);
);
Write( "\!n",good," of ", N Items( condat )," names resolved." );Le bit suivant concerne les informations d'identification : le ppk est la clé privée de la machine JMP. La partie publique de cette clé se trouve sur chaque machine du cluster. Plink a besoin de ce bit privé pour signer la connexion SSH à la machine distante. Le code if(0) est là pour que je puisse expérimenter sur toutes les machines. À la fin de la journée, je pourrais leur envoyer à tous la commande d'arrêt. Soyez prudent avec la commande shutdown si la liste condat contient la machine Linux qui exécute la VirtualBox qui exécute ce script JSL (ligne 9 ci-dessus.)
remoteuser = "c@";
localkey = " -i c:\Users\v1\Desktop\keys\private.ppk ";
// test code examples...
if( 0, // sequential test using command
//command = " \!"python3.6 -m pip install numpy";
//command = " \!"python3.6 -m pip --version";
//command = " \!"python3.6 -m pip install -U pip";
//command = " \!"python3.6 --version";
command = " \!"hostname -A";
command = " \!"top -b -n 2|head -n 15";
// check if the common nfs mount is set up
//command = " \!"ls /mnt/nfsc"; // fstab: closet.local:/mnt/nfs /mnt/nfsc/ nfs defaults 0 0
// cycle the power; -r means restart. for this to work, must be sudo and not need password
//command = " \!"sudo shutdown -r now"; // visudo: c ALL=(ALL) NOPASSWD: /sbin/poweroff, /sbin/reboot, /sbin/shutdown
//command = " \!"sudo shutdown now"; // visudo: c ALL=(ALL) NOPASSWD: /sbin/poweroff, /sbin/reboot, /sbin/shutdown
// info about number of cores
//command = " \!"lscpu | egrep '^(Model name|Thread|CPU\(s\):)'";
// info about memory
//command = " \!"free -h";
for( iCondat = 1, iCondat <= N Items( condat ), iCondat += 1,
result = runprogram(
executable( "plink" ),
options(
localkey ||
/**/ " -hostkey " || condat[iCondat][hostkey] || /* this line hides the host key */
" " || remoteuser || condat[iCondat][server] || command
),
write function( Function( {this}, this << Write( "n\r\n" ) ) ),
readfunction( "blob" )
);
result = Regex( Blob To Char( result ), "(\r\n)+", "\r\n", GLOBALREPLACE );
Write( "\!n\!n", iCondat, " ", condat[iCondat][server], " ", result );
Wait( 0 );
);
Write( "\!n\!ntest done" );
);Le bit suivant concerne l'exécution simultanée de plusieurs copies de RunProgram. ImgNum est la clé de chacun de ces tableaux associatifs et correspond au numéro de l'image sur laquelle le programme Python doit travailler. ImgNumToRP contient un objet RunProgram en cours d'exécution. Cet objet est déplacé vers imgNumToRPdone lorsqu'il se termine. imgNumToResult contient la sortie texte du programme Python, imgNumToMachine l'identifiant de la machine et imgNumToLabel une étiquette pour l'image. (Les étiquettes sont appliquées plus tard, plus ci-dessous.)
imgNumToRP = [=> ]; // holder for active RunProgram objects
imgNumToRPdone = [=> ]; // holder for RunProgram objects that need a moment to terminate
imgNumToRPresult = [=> ""]; // result strings from RunProgram objects
imgNumToMachine = [=> 0]; // which machine is this running on?
imgNumToLabel = [=>""]; // f(x) labels for graphs
nextImageNumberToMakePicture = 0;
// also: keep track of zeros from previous image
La fonction unstick gère les RunPrograms qui n'appellent pas la ReadFunction une dernière fois après la fin de RunProgram.
unstick = function({},{imagenumber},
for( imagenumber = imgNumToRP << first, !Is Empty( imagenumber ), imagenumber = imgNumToRPnext,
imgNumToRPnext = imgNumToRP << next( imagenumber );
if( (imgNumToRP[imagenumber]) << canread == 0 & (imgNumToRP[imagenumber]) << canwrite == 0 & (imgNumToRP[imagenumber]) << isreadeof == 1,
imgNumToRPresult[imagenumber] = "\!n" || char(imagenumber) || ": " || imgNumToRPresult[imagenumber] || "\!n"; // NoInterleave save this with count and basin
imgNumToRPdone[imagenumber] = imgNumToRP[imagenumber];
imgNumToRP << removeitem( imagenumber );
//if(imgNumToMachine[imagenumber]!=machine,throw("b machine not match"||char(machine)||","||char(imgNumToMachine[imagenumber])));
condat[imgNumToMachine[imagenumber]][running] = condat[imgNumToMachine[imagenumber]][running] - 1;
imgNumToMachine[imagenumber] = 0;
try(
write("\!nunsticking ",imagenumber);
process( imagenumber ); // must be after...
,
show(exception_msg);
);
Break(); // restart with <<FIRST></FIRST>La fonction findmachine ci-dessous recherche une machine qui n'a pas démarré trop de sessions SSH. La fonction unstick() est appelée périodiquement pour rechercher les 10 % de cas bloqués.
findMachine = function( {},{iCondat,n}, // find a machine that is not filled up
n=0;
while( 1,
for( iCondat = 1, iCondat <= N Items( condat ), iCondat += 1,
if( condat[iCondat][running] < condat[iCondat][cpus],
Return( iCondat )
)
);
if(mod(n+=1,60)==0,
write("wait for machine...");
// look for stuck
unstick();
);
Wait( 1 ); // a full second is fine; the processes run for 10s of seconds
)
);La fonction de processus est appelée lorsque JMP détecte qu'un RunProgram est terminé. Les programmes Python écrivent "sentinelle" à la fin de leur sortie et ce code le vérifie comme une vérification brute que tout fonctionne correctement. Bindata est le fichier binaire que les programmes Python produisent dans le répertoire de la ferme sur le disque NFS. Windows et JMP voient le disque NFS au chemin \\VBOXSVR mais Linux et les programmes Python le voient au chemin //mnt. Ce code allait produire les images à partir des données binaires mais cela a été séparé en une étape de post-traitement. Maintenant, il rapporte simplement le nombre de zéros découverts pour la formule. Le format binaire est expliqué dans le code python ; JMP le lit et le charge avec blobtomatrix, utilisé deux fois ici. La relance à la fin fonctionne probablement si elle n'est pas commentée et si nécessaire ; J'ai oublié de le décommenter et je n'en ai pas eu besoin.
// process output from one thread in the farm
process = function( {imagenumber},
{bindata, start, stop, balanceBytes, matdata, countArray,basinArray,zeros},
//show("process",imagenumber);
if( Contains( imgNumToRPresult[imagenumber], "sentinel" ),
// process a file. maybe not this one if it is out of order... FIXME...nextImageNumberToMakePicture
bindata = Load Text File( "\\VBOXSVR\nfsc\farm\binary" || char(imagenumber) || ".generated", blob() );
matdata = blobtomatrix(
bindata,
"uint",
2,
"little",
2
);
countArray = J( yresolution, xresolution, 0 );
basinArray = J( yresolution, xresolution, 0 );
countArray[0] = matdata[1::nrows(countArray)*ncols(countArray), 1];
basinArray[0] = matdata[1::nrows(countArray)*ncols(countArray), 2];
balanceBytes = length(bindata) - 4*nrows(countArray)*ncols(countArray);
zeros = blobtomatrix(blobpeek(bindata,length(bindata)-balanceBytes),"float",8,"little",4);
write("\!nimagenumber=",imagenumber, if(!ismatrix(zeros)," zeros ***"," zeros="||char(nrows(zeros))));
deletefile("\\VBOXSVR\nfsc\farm\testpython"||char(imagenumber)||".py"); // source no longer needed
sentinels += 1;//
, // else
Write( imgNumToRPresult[imagenumber] );
Write( "\!nrelaunch disabled" || char(imagenumber) );
//launch( imagenumber );
)
);La fonction de lancement démarre un autre RunProgram. Il utilise findmachine pour attendre une machine disponible et fournit un fichier .py personnalisé et des paramètres personnalisés à la commande Python. Le fichier .py est créé plus bas. Le lancer à la ligne 5 ci-dessous est un test de santé mentale ; le JSL est assez compliqué, j'ai senti que j'avais besoin d'écrire du code supplémentaire pour surveiller les erreurs de logique.
// launch a thread in the farm
launch = function( {imagenumber},{command,machine},
machine = findMachine();
//show("launch",imagenumber,machine);
if(imgNumToMachine[imagenumber]!=0,throw("image number already on machine"||char(imagenumber)||","||char(machine)));
imgNumToMachine[imagenumber] = machine;
condat[machine][running] = condat[machine][running] + 1;
//command = "\[python3 /mnt/nfsc/farm/testpython]\"||char(imagenumber)||"\[.py]\";
command = "\[hostname -A;python3 /mnt/nfsc/farm/testpython]\"||char(imagenumber)||"\[.py]\";
imgNumToRP[imagenumber] = runprogram(
executable( "plink" ),
options(
localkey || /* -v */ " -hostkey " || condat[machine][hostkey] || " " || remoteuser || condat[machine][server] || " " || command || //
" -o " || Char( 0 ) || // offset
" -s " || Char( 1 ) || // step
" -x " || Char( xresolution ) || //
" -y " || Char( yresolution ) || //
" -t " || Char( ytop ) || //
" -l " || Char( xleft ) || //
" -r " || Char( xright ) || //
" -b " || Char( ybottom ) || //
" -f " || "/mnt/nfsc/farm/binary" || char(imagenumber) || ".generated" || //
""
),Le mot clé de paramètre rend son argument disponible à la fonction read. C'est ainsi que vous pouvez avoir plusieurs RunPrograms en parallèle et que les fonctions de lecture savent de quoi il s'agit.
Parameter( evallist({imagenumber,machine}) ),La fonction de lecture capture la sortie et essaie (environ 90 % du temps avec succès) de détecter que RunProgram est terminé. Étant donné que RunProgram se termine toujours et s'exécute toujours lorsque ce code s'exécute, le handle est déplacé vers imgNumToRPDone au lieu d'être simplement supprimé. Cela lui permet de finir proprement. Parm, ci-dessous, reçoit la valeur du paramètre, ci-dessus. Il y a un autre contrôle de santé mentale (lancer) ci-dessous.
readfunction(
function( {this, parm}, {txt, imagenumber,machine},
{imagenumber,machine} = parm ;
//show("readfunc",imagenumber,machine);
txt = "";
//write("readfun("||char(imagenumber)||")... ");
While( this << canread, txt = txt || (this << read) );
if( Length( txt ) > 0,
//show(txt);
imgNumToRPresult[imagenumber] = imgNumToRPresult[imagenumber] || txt
);
if( this << isreadeof & imgNumToRP << Contains( imagenumber ), // this *almost* always works, about 1:100 fail to get the readeof in the readrunction
imgNumToRPresult[imagenumber] = "\!n" || char(imagenumber) || ": " || imgNumToRPresult[imagenumber] || "\!n"; // NoInterleave save this with count and basin
process( imagenumber );
imgNumToRPdone[imagenumber] = imgNumToRP[imagenumber];
imgNumToRP << removeitem( imagenumber );
if(imgNumToMachine[imagenumber]!=machine,throw("a machine not match"||char(machine)||","||char(imgNumToMachine[imagenumber])));
condat[machine][running] = condat[machine][running] - 1;
imgNumToMachine[imagenumber] = 0;
);
)
)
);
(imgNumToRP[imagenumber]) << writeeof;
);La fonction makePythonCode génère le fichier .py avec l'équation oneStep intégrée. Le programme python inclut une fonction init() qui analyse la ligne de commande pour obtenir des informations sur la taille de l'image et le nom du fichier de sortie. Les données sont transmises de deux manières : paramètres de ligne de commande et intégrées dans des équations. Peut-être que toutes les informations auraient pu être intégrées, mais s'il n'y a aucune raison (comme compiler l'équation oneStep dans le code) d'intégrer des paramètres dans le code, vous n'auriez besoin que d'une copie du code pour tous les RunPrograms s'il était entièrement contrôlé par les paramètres de ligne de commande. Le code Python est écrit de manière à résoudre une partie du problème si je voulais que plusieurs ordinateurs travaillent sur une seule image. Il est globalement plus rapide de créer une image complète en une seule fois et de ne pas recombiner les parties ultérieurement.
makePythonCode = function({imagenumber,oneStep},{filename,python},
python =
"\[
import cmath, math, sys, getopt
import numpy as np
def main():
offset, step, ixdim, iydim, top, bottom, left, right, filename = init() # NoInterleave offset will be 0 and step will be 1
iymax = iydim - 1
ixmax = ixdim - 1
yrange = top - bottom
xrange = right - left
result = np.zeros( (int(math.ceil((iydim*ixdim-offset)/step)),2),dtype=np.int16) # count, basin # NoInterleave simpler, no change required though
iresult = 0
ix = offset # NoInterleave like iy, now begins at zero
iy = 0
while ix > ixmax: # normalize before loop starts
ix = ix - ixdim
iy = iy + 1
y = bottom + yrange * iy / iymax # after iy is normalized, before loop
roundedComplexToBasinId = {}
basinIdToBasinSize = {} # NoInterleave this might be useful, but needs to be appended to file (now a constant size for main chunk)
nextBasinId = 1 # reserve basin 1 for "did not converge"
while iy <= iymax:
x = left + xrange * ix / ixmax
xOld = complex(x,y) #xOld - (xOld ** 5 - 1) / (5 * xOld ** 4)
count = 0
try: ## catch zero divide, rarely
for irun in range(1,1000): # NoInterleave prev stopped at 100
###xNew = xOld - (xOld ** 5 - 1) / (5 * xOld ** 4)
]\" || // maintain the indent and newlines!La ligne suivante est le bit intégré… OneStep est une chaîne de caractères JSL d'une instruction Python valide qui définit xNew à partir de la valeur xOld et de l'équation de Newton, comme le commentaire ### deux lignes ci-dessus.
oneStep || "\[
count += 1
delta = abs(xNew - xOld)
dist = abs(xOld)
if delta < 1e-15 * dist:
break
xOld = xNew
xOld = complex( round(xOld.real,10), round(xOld.imag,10) )
if irun >= 999:
basinId = 1 # arbitrary, but available always. did not converge.
else:
try: ## catch new basin discovered, rarely
basinId = roundedComplexToBasinId[xOld]
except KeyError as e:
nextBasinId += 1
roundedComplexToBasinId[xOld] = nextBasinId
basinId = nextBasinId
basinIdToBasinSize[basinId] = 0
### print(offset, count, xOld, basinId)
except ZeroDivisionError as e:
basinId = -1
### print(offset, e)
try:
result[iresult,0] = count
except:
print("oops",ix,iy)
if basinId >= 0:
try:
basinIdToBasinSize[basinId] += 1
except KeyError as e:
basinIdToBasinSize[basinId] = 1 # happens with basinId==1
result[iresult,1] = basinId
iresult += 1
ix = ix + step # advance
### offset = offset + step # needed for return value
while ix > ixmax: # re-normalize
ix = ix - ixdim
iy = iy + 1
y = bottom + yrange * iy / iymaxVoici le code python qui écrit un fichier binaire d'entiers et de doubles. Les doubles suivent les entiers pour permettre à JMP de saisir facilement tout avec blobToMatrix.
# output results
fileobj = open(filename, mode='wb') # NoInterleave maybe prepend a descriptor
result.tofile(fileobj) # NoInterleave append after this constant size block
# NoInterleave append roundedComplexToBasinId and basinIdToBasinSize rather than printing them, they should be same dimension
# 4 8-byte doubles: real,imag,basin,size
temp = np.zeros( (len(roundedComplexToBasinId),4),dtype=np.float64)
itemp = 0
for iz in roundedComplexToBasinId.keys():
temp[itemp,0] = iz.real
temp[itemp,1] = iz.imag
temp[itemp,2] = roundedComplexToBasinId[iz]
temp[itemp,3] = basinIdToBasinSize[roundedComplexToBasinId[iz]]
itemp += 1
temp.tofile(fileobj)
fileobj.close()
print("sentinel"); # so JSL can test if this ran to the end...
def init():
offset = 0
step = 1
xdim = 3
ydim = 2
top = -1.0
bottom = 1.0
left = -1.0
right = 1.0
filename = "temp.binary"
options, remainder = getopt.getopt(sys.argv[1:], 's:x:y:o:t:b:l:r:f:')
for opt, arg in options:
if opt == '-o':
offset = int(arg) # offset is the starting point
elif opt == '-s':
step = int(arg) # step to the next point
elif opt == '-x':
xdim = int(arg) # xdim field int size
elif opt == '-y':
ydim = int(arg) # ydim field int size
elif opt == '-t':
top = float(arg) # top of field user coords
elif opt == '-b':
bottom = float(arg) # bot field user coords
elif opt == '-l':
left = float(arg) # left field user coords
elif opt == '-r':
right = float(arg) # right field user coords
elif opt == '-f':
filename = arg # write binary file
print(" xdim=",xdim,
" ydim=",ydim,
" offset=",offset,
" step=",step,
" top=",top,
" bottom=",bottom,
" left=",left,
" right=",right)
return offset, step, xdim, ydim, top, bottom, left, right, filename
if __name__ == "__main__":
main()
]\";Le programme Python est enregistré dans le répertoire NFS (comme Windows le connaît) avec un numéro unique. Le nom de fichier renvoyé n'est pas utilisé ; la version Linux du nom de fichier est créée par la fonction de lancement et transmise via Plink. (Nous sommes de retour dans JSL maintenant...)
filename = Save Text File( "\\VBOXSVR\nfsc\farm\testpython"||char(imagenumber)||".py", python ); // NoInterleave - unique file name for each picture because the equation changes
return(filename);
);
La boucle ci-dessous a généré les 10 000 images.
for( imagenumber = 0, imagenumber < 10000, imagenumber += 1, //<<<<< frame loop <<<<<<
write("\!nloop",imagenumber);
// for X^5 - 1 there are 5 roots; 1 real and 4 complex
// the python expression for x(n+1) = x(n) - F / F' is
// xNew = xOld - (xOld ** 5 - 1) / (5 * xOld ** 4)La formule et la conversion de OneStep en une instruction Python - cela semble compliqué car il remplace les valeurs réelles dans les paramètres F1...F6 et E1...E6. Caché à la fin est le changement de ^ de JSL à ** de Python, qui se fait dans le texte, pas dans l'expression, car JSL n'a aucun moyen de représenter ** comme une expression. Les F1...E6 sont traités comme des expressions. Sin( imagenumber / 200 ) correspond à environ 10 ondes sinusoïdales complètes de -1 à 1. Les ondes sinusoïdales ont des fréquences différentes en raison des différents diviseurs. Les facteurs se déplacent en douceur de -2 à 2 et les exposants se déplacent en douceur autour de leur cible 1...6 de +/- 1. Ainsi, l'équation d'une image à l'autre est assez similaire, mais varie un peu de manière imprévisible.
FF = Derivative( Name Expr( F ), xOld );
// onestep will be inserted into python below
oneStep = substitute( // "xNew = xOld - (xOld ** 5 - 1) / (5 * xOld ** 4)"
char( // "xNew = xOld - (xOld ^ 5 - 1) / (5 * xOld ^ 4)"
substitute( // char() adds the parens around F and FF
Expr( xNew = xOld - F / FF ), // build the statement with symbols, then
Expr( F ), Name Expr( F ), // replace F and FF symbols with the content
Expr( FF ), Name Expr( FF ), // of the function and derivative
Expr( F6 ), round(2 * Sin( imagenumber / (200+5) ) / 2,3),
Expr( F5 ), round(2 * Sin( imagenumber / (200+7) ) / 2,3),
Expr( F4 ), round(2 * Sin( imagenumber / (200+11) ) / 2,3),
Expr( F3 ), round(2 * Sin( imagenumber / (200+13) ) / 2,3),
Expr( F2 ), round(2 * Sin( imagenumber / (200+17) ) / 2,3),
Expr( F1 ), round(2 * Sin( imagenumber / (200+19) ) / 2,3),
Expr( E6 ), round(6.0 + Sin( imagenumber / (200+23) ) / 2,3),
Expr( E5 ), round(5.0 + Sin( imagenumber / (200+29) ) / 2,3),
Expr( E4 ), round(4.0 + Sin( imagenumber / (200+31) ) / 2,3),
Expr( E3 ), round(3.0 + Sin( imagenumber / (200+37) ) / 2,3),
Expr( E2 ), round(2.0 + Sin( imagenumber / (200+41) ) / 2,3),
Expr( E1 ), round(1.0 + Sin( imagenumber / (200+47) ) / 2,3)
) // sub F,FF
), // char
"^", "**"
); // substitute ^ in JSL to ** in pythonLe tableau associatif imgNumToLabel peut être utilisé ultérieurement pour étiqueter les images. Simplifier xAncien->X...
imgNumToLabel[imagenumber]=substitute(regex(oneStep,"xNew = xOld - \((.*?)\) / ","\1")," ** ","^"," * ","·","xOld","X");
write("\!n",imgNumToLabel[imagenumber]);
//"/mnt/nfsc/farm/binary" || char(imagenumber) || ".generated"Le test d'un fichier existant a été utilisé pendant le développement pour accélérer les choses lorsque le programme devait être redémarré. Cela peut ne pas fonctionner correctement, il n'a pas été utilisé pour le dernier run (trois jours !)
if(fileexists("\\VBOXSVR\nfsc\farm\binary" || char(imagenumber) || ".generated"),
write("Exists ","\\VBOXSVR\nfsc\farm\binary" || char(imagenumber) || ".generated");
continue();
);
makePythonCode(imagenumber,oneStep);Cela ressemble à du code abandonné déplacé.
n = 0; // eventually this becomes the offset // NoInterleave offset is always 0
sentinels = 0; // NoInterleave "sentinel" might be all it returns
//write("Launching... ");C'est là que commencent les 10 000 RunPrograms. 50 ou 60 fonctionneront en même temps en utilisant les valeurs de condat, bien au-dessus.
launch( imagenumber );
write("Launched. ")
);La conclusion... attendez que les derniers RunPrograms se terminent, unstick peut en trouver un ou deux.
while( N Items( imgNumToRP ) > 0,
Wait( .1 );
unstick();
);
write("\!ndone");Les étiquettes sont enregistrées pour le générateur d'images, plus tard. JSL/Python ci-dessus a déjà enregistré des données binaires qui doivent encore être converties au format .png.
savetextfile("$desktop/imgNumToLabel.txt",char(imgNumToLabel));Le Rpdone contient 10 000 RunPrograms terminés. Laisser aller.
imgNumToRPdone = 0;// clean out // NoInterleave picture generation happens after all count and basin are created
Le générateur d'images peut être exécuté maintenant ou plus tard...
include("pythonFarmDisplay.jsl");C'est la fin du code Beowulf. Le générateur d'images aurait pu l'utiliser, mais ne l'a pas fait. Je veux que la vidéo utilise les mêmes couleurs d'une image à l'autre. Les couleurs représentent les bassins d'attraction, le zéro vers lequel tout point du bassin est attiré doit avoir la même couleur. Mais l'équation change un peu à chaque image, donc les zéros se déplacent. Le générateur d'images essaie de faire correspondre les zéros de l'image précédente avec l'image suivante. Parfois, ils sautent, mais c'est assez bon.
Le JSL pour faire les photos suit. Faites d'abord une table de couleurs pour les bassins. Ce sont la version la plus claire possible de la couleur; ils sont assombris plus tard par le grand nombre d'itérations.
colors=
hlscolor(0,L,1)|| // dark red
hlscolor(.05,L,0.8)|| // red
hlscolor(.15,L,0.9)|| // orn
hlscolor(.35,L,1.0)|| // yel
hlscolor(.50,L,1.0)|| // grn
hlscolor(.65,L,1.0)|| // cyn
hlscolor(.80,L,0.7)|| // blu
hlscolor(.95,L,1.0)|| // mag
hlscolor(1.0,L,0.0); // gry
redmap=j(ncols(colors),1,0);
greenmap=j(ncols(colors),1,0);
bluemap=j(ncols(colors),1,0);
for(i=1,i<=ncols(colors),i+=1,
{r,g,b}=colortorgb(colors[i]);
m = max(r,g,b);
redmap[i] = r/m;
greenmap[i] = g/m;
bluemap[i] = b/m;
);Les étiquettes d'image enregistrées précédemment doivent être rechargées, les dimensions de l'image doivent correspondre. Ce code doit être le même que précédemment et doit avoir été extrait dans un fichier de paramètres commun. Trop tard maintenant!
imgNumToLabel=Include( "$desktop/pythonFarmDisplayLabels.jsl" );
scale = 2 / 1;
xresolution = 1920 * scale;
yresolution = 1080 * scale;
// keep aspect ratio by chopping shorter dimension...
xcenter = -1;
ycenter = 0;
delta = 3;
ytop = ycenter - delta * yresolution / xresolution;
xleft = xcenter - delta;
xright = xcenter + delta;
ybottom = ycenter + delta * yresolution / xresolution;
Les couleurs des bassins sont censées suivre des images/équations légèrement différentes. C'est là que vivent les zéros de l'image précédente.
prevzeros = []; // what came before the first one
Le disque NFS était bien trop petit. Il y a 309 Go d'informations sur les pixels non compressés ici.
//filepat = "\\VBOXSVR\nfsc\farm\binary^imagenumber^.generated";
// manually moved the files (running out of space)...
filepat = "f:\newtonmethod\generated\binary^imagenumber^.generated";
countarray = J( yresolution, xresolution, 0 );
basinarray = J( yresolution, xresolution, 0 );
Les fichiers binaires seront décompressés en deux matrices JSL : un décompte du nombre d'itérations nécessaires au pixel et un petit entier identifiant le bassin du pixel. Il y a un peu plus de données binaires après cela (les octets d'équilibre) qui amènent le numéro du bassin à une coordonnée complexe ; la coordonnée complexe est ce qui correspond d'une image à l'autre.
For( imagenumber = 1, imagenumber < 10000, imagenumber += 1,
file = Eval Insert( filepat );
If( File Exists( file ),
bindata = Load Text File( file, blob() );
matdata = Blob To Matrix( bindata, "uint", 2, "little", 2 );
countArray[0] = matdata[1 :: yresolution * xresolution, 1];
basinArray[0] = matdata[1 :: yresolution * xresolution, 2];
balanceBytes = Length( bindata ) - 4 * yresolution * xresolution;
zeros = Blob To Matrix( Blob Peek( bindata, Length( bindata ) - balanceBytes ), "float", 8, "little", 4 );
//Show( zeros );
if(zeros[1, 3]==2,
zeros = [9e99 9e99 1 0] |/ zeros; // python code should have put this in
);
// re-assign the basin numbers to line up with prevzeros
If( N Rows( prevzeros ) == 0,
prevzeros = zeros
);
Le KDTable trouve les voisins les plus proches ; le plan complexe est un beau voisinage.
kdt = KDTable( prevzeros[0, 1 :: 2] ); // 1::2 means the real and imag cols
newbasin = (1 :: N Rows( zeros )) + 10000;
toomany = N Rows( zeros ) > N Rows( prevzeros );
For( i = 1, i <= N Rows( zeros ), i += 1,
If( zeros[i, 3] != i,
Throw( "basin index assumption wrong" )
); // zeros are sorted by basin number
{indx, dist} = kdt << knearestrows( 2, zeros[i, 1 :: 2] );
If( N Cols( dist ) > 1 & dist[1] < dist[2] / 10, // assume 10X closer is a good match
newbasin[i] = indx[1];
kdt << removerows( indx[1] );//
, /*else*/ !toomany, // if there are toomany new zeros, don't assign them this way
newbasin[i] = indx[1];
kdt << removerows( indx[1] );
);
);<
/*
zerobasin newbasin
1 4 change 1s to 4s in basinArray and zeros
2 3 2s to 3s
3 10001 too many case! no basin assigned. becomes 7, below
4 2
5 5
6 6
7 1
*/
unassigned = Loc( newbasin > 10000 );
If( N Rows( unassigned ) > 0,La nouvelle image a plus de zéros que l'ancienne image, placez-les dans le dernier bassin.
newbasin[unassigned] = newbasin[unassigned] - 10000 + N Rows( zeros )
);
basinArray = Shape( newbasin[basinArray], N Rows( basinArray ), N Cols( basinArray ) );
zeros[0, 3] = Shape( newbasin[zeros[0, 3]], N Rows( zeros ), 1 );
ri = Rank Index( zeros[0, 3] );
zeros[0, 0] = zeros[ri, 0]; // zeros are sorted by basin number
Les zéros sont triés au-dessus ; la ligne 1 obtient la couleur 1, etc.
prevzeros = zeros;
Les images avec un nombre élevé sont supprimées pour rendre la vidéo plus fluide.
countarray[loc(countarray>200)] = 200; // 150 is typical big count, 999 is not converged
rangecount = Range( countarray );
countarray=countarray-3; // 3 or 4 is typical minimum count, make it bright (0)
countarray[loc(countarray<0)]=0;
Le GrayArray est mis à l'échelle de l'obscurité à partir des comptes
grayArray = (1 - (countArray /*- rangeCount[1]*/) / (200 /*rangeCount[2] - rangeCount[1]*/)) ^ 8;
rangebasin = range(basinarray);
//basinArray = basinArray / (7);
nc=ncols(colors);
basinarray[loc(basinarray>nc)]=ncols(colors); // if these exist they could be mapped to gray, perhaps
basinarray[loc(basinarray<1)]=ncols(colors); // ditto
Le gray * redmap[basinarray] (etc) est un composant de couleur RVB mis à l'échelle
red = grayArray :*</img> shape(redmap[basinarray], N Rows( basinArray ), N Cols( basinArray ) );//(1 + Sin( 2 * Pi() * basinArray + (0 * 2 * Pi() / 3) )) / 2;
green = grayArray :*</img> shape(greenmap[basinarray], N Rows( basinArray ), N Cols( basinArray ) );//(1 + Sin( 2 * Pi() * basinArray + (1 * 2 * Pi() / 3) )) / 2;
blue = grayArray :*</img> shape(bluemap[basinarray], N Rows( basinArray ), N Cols( basinArray ) );//(1 + Sin( 2 * Pi() * basinArray + (2 * 2 * Pi() / 3) )) / 2;
Les correctifs pour le dépassement/sous-dépassement des valeurs de couleur 0..1…
red[Loc( red < 0 )] = 0;
red[Loc( red > 1 )] = 1;
green[Loc( green < 0 )] = 0;
green[Loc( green > 1 )] = 1;
blue[Loc( blue < 0 )] = 0;
blue[Loc( blue > 1 )] = 1;
Les marqueurs + sont ajoutés aux emplacements zéro… brut, mais fonctionne. Ne dessinez pas sur le bord, l'indexation matricielle fait des lignes ou des colonnes entières pour un indice de zéro !
crosssize = floor(yresolution / 100);
For( iz = 1, iz <= N Rows( zeros ), iz += 1,
real = Round( xresolution * Interpolate( Num( zeros[iz, 1] ), xleft, 0, xright, 1 ) );
imag = Round( yresolution * Interpolate( Num( zeros[iz, 2] ), ytop, 0, ybottom, 1 ) );
if(crosssize+5 < imag < nrows(red)-crosssize-5 & crosssize+5 < real < ncols(red)-crosssize-5,
red[imag - 1, real - crosssize :: real + crosssize] = 0.5;
red[imag + 1, real - crosssize :: real + crosssize] = 0.5;
red[imag - 2, real - crosssize :: real + crosssize] = 0.5;
red[imag + 2, real - crosssize :: real + crosssize] = 0.5;
red[imag - crosssize :: imag + crosssize, real - 1] = 0.5;
red[imag - crosssize :: imag + crosssize, real + 1] = 0.5;
red[imag - crosssize :: imag + crosssize, real - 2] = 0.5;
red[imag - crosssize :: imag + crosssize, real + 2] = 0.5;
green[imag - 1, real - crosssize :: real + crosssize] = 0.5;
green[imag + 1, real - crosssize :: real + crosssize] = 0.5;
green[imag - 2, real - crosssize :: real + crosssize] = 0.5;
green[imag + 2, real - crosssize :: real + crosssize] = 0.5;
green[imag - crosssize :: imag + crosssize, real - 1] = 0.5;
green[imag - crosssize :: imag + crosssize, real + 1] = 0.5;
green[imag - crosssize :: imag + crosssize, real - 2] = 0.5;
green[imag - crosssize :: imag + crosssize, real + 2] = 0.5;
blue[imag - 1, real - crosssize :: real + crosssize] = 0.5;
blue[imag + 1, real - crosssize :: real + crosssize] = 0.5;
blue[imag - 2, real - crosssize :: real + crosssize] = 0.5;
blue[imag + 2, real - crosssize :: real + crosssize] = 0.5;
blue[imag - crosssize :: imag + crosssize, real - 1] = 0.5;
blue[imag - crosssize :: imag + crosssize, real + 1] = 0.5;
blue[imag - crosssize :: imag + crosssize, real - 2] = 0.5;
blue[imag - crosssize :: imag + crosssize, real + 2] = 0.5;
red[imag, real - crosssize :: real + crosssize] = 1;
red[imag - crosssize :: imag + crosssize, real] = 1;
green[imag, real - crosssize :: real + crosssize] = 1;
green[imag - crosssize :: imag + crosssize, real] = 1;
blue[imag, real - crosssize :: real + crosssize] = 1;
blue[imag - crosssize :: imag + crosssize, real] = 1;
);
);
label = imgNumToLabel[imagenumber];L'étiquette est ajoutée en bas. Faites un petit bitmap de l'étiquette, puis copiez-le dans l'image plus grande.
if(1,
//label = char(imagenumber)||" zeros="||char(N Rows( zeros ))||" "||label;
labelpic = textbox(label,setwrap(3000),<<SETFONTSIZE></SETFONTSIZE>
L'image est faite ici. Décidez ensuite si elle sera conservée dans la vidéo en la nommant.
img = New Image( "rgb", {red, green, blue} );
if(6<=nrows(zeros)<=8 & 2<=rangebasin[1]<=4 & 6<=rangebasin[2]<=8 & 3<=rangecount[1]<=4 & rangecount[2]<200,
img << saveimage( "f:/newtonmethod/png3/keeper" || Char( 1e6 + imagenumber ) || ".png", "png" );
,
img << saveimage( "f:/newtonmethod/png3/zapped" || Char( 1e6 + imagenumber ) || ".png", "png" );
);
//img<<SCALE></SCALE>C'est ça! Voici la vidéo. J'ai décidé de ne faire que la version scientifique, avec des cadres sombres scintillants là où ils se produisent et une légende occupée ci-dessous. Une jolie version pourrait être plus douce et encore plus hypnotique. C'est une vidéo 4K, si vous avez le matériel et la vitesse de connexion, sinon youtube la réduira pour vous.
Il reste deux heures de 2020 au moment où je poste ceci. J'ai vraiment hâte d'être en 2021.