This website uses Cookies. Click Accept to agree to our website's cookie use as described in our Privacy Policy. Click Preferences to customize your cookie settings.
- JMP User Community
- :
- Blogs
- :
- Uncharted
- :
- Multiple precision arithmetic with Python
Turn on suggestions
Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.
Article Options
- Subscribe to RSS Feed
- Mark as New
- Mark as Read
- Bookmark
- Subscribe
- Printer Friendly Page
- Report Inappropriate Content
Super User
Multiple precision arithmetic with Python
Created:
Jan 1, 2022 06:59 PM
Playing with JMP + Python + gmpy2 (python library for multiple precision arithmetic) + mpfr (multiple precision real) + mpc (multiple precision complex). Using the same Beowulf cluster descibed in Beowulf, Newton, and Mr Hanson . It is enormously slow, even with 8 machines/~44 CPUs available. In the video you can see a window monitoring 100 threads in a 10x10 grid; each thread is computing the pixels for a small patch of a zoomed in Mandelbrot. With normal (fast) arithmetic, getting past a zoom of 1e-13 or so is impossible because there are only about 15 significant digits. In the video, at the end, the zoom is around 1e-21. It is getting really slow because the iteration count is 200,000 for the hardest pixels, and probably close to that for most of them. Thanks to online clock for the real time clock in the bottom right.
Here's the JSL I'm playing with; it is set up for my environment and can't possibly "just work" but it might be useful for an example of
- many RunPrograms working in parallel
- using a cluster of Linux machines running Python
- using the Python Multiple Precision library (you might not need that, apparently Python natively supports really big integers.)
- putting thumbnails in a data table (better places to see that than here maybe)
The multiple precision numbers are represented as strings in JMP. The Python program has a tiny calculator JMP uses when zooming in; JMP tells the calculator the 0.0 to 1.0 coordinates of a selection rectangle in the graph box (no axes visible, but still a graph box) and Python does the mpfr calculation and returns the long string of resulting digits.
// 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
condat = { // load balence with audacity running on thinkb and using 100...good enough...
// the 3,3,4,4,4,5,6,8 balance is not just about machine speed, but about the
// poor algorithm for asigning machines that uses up line 1, then 2, ... then 8
// then sets back to line 1. for 100 that means
// 3+3+4+4+4+5+6+8+ 3+3+4+4+4+5+6+8 +3+3+4+4+4+5+3
// 9 9 12 12 12 15(b) 17(c) 19(d)
// then 'think' machines are a little faster than the 9010/7010 machines
// and twice as fast as the 790 machines. (not 8/3, but maybe 17/9 or 19/9)
{"4e:c3:8a:3b:e7:d4:bd:0d:85:62:36:40:17:15:b4:a3", "dell790a.local", 3, 0, 22},
{"00:c2:37:09:dc:d7:cf:fb:33:e3:14:e3:d6:fb:39:be", "dell790b.local", 3, 0, 14},
{"54:c8:1a:6a:50:e3:51:02:ed:1b:f8:fa:bd:f4:e4:17", "dell7010a.local", 4, 0, 30},
{"66:0d:ca:3c:1f:11:eb:45:d4:92:d0:0b:ae:01:b8:6b", "dell7010b.local", 4, 0, 30},
{"94:84:2f:a5:56:2d:db:16:b5:80:a9:0d:87:07:ee:53", "dell9010.local", 5, 0, 30},
{"b2:b7:b3:c3:83:63:11:e7:56:e7:c9:c2:e6:fc:f0:c9", "thinkb.local", 6, 0, 30},
{"2b:ca:fc:45:f0:64:cd:ac:70:ce:93:73:e8:84:2d:75", "thinkc.local", 6, 0, 30},
{"62:b7:67:23:8d:86:fe:80:1a:62:ec:d8:41:1b:b4:1b", "thinkd.local", 6, 0, 30}
//{"9a:e1:a7:05:ed:f7:d0:4c:1d:c6:cb:d2:db:f3:02:c4", "d790main.local", 0, 0, 30}, // host for vm
//{"14:26:bc:41:48:ef:43:31:61:29:e5:a6:3f:d5:43:07", "closet.local", 0, 0, 14} // tv capture
};
machinename = {};
for(i=1,i<=nitems(condat),i+=1,machinename[i]=regex(condat[i][server],"[^.]+"));
good = 0;
For( iCondat = 1, iCondat <= N Items( condat ), iCondat += 1,
condat[iCondat][running] = 0;
ai = Get Addr Info( condat[iCondat][server] );
For( j = 3, j <= N Items( ai ), j += 1,
If( "PF_INET" == ai[j][1], // needs [ ] wrapped around "PF_INET6" addrs
Try(
ip = ai[j][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] ) )
)
)
);
Wait( .01 );
);
If( good != N Items( condat ),
Throw( "did not find all condat" );
Wait( .1 );
);
Write( "\!n", good, " of ", N Items( condat ), " names resolved." );
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 -haltr /mnt/nfsc/hi"; // fstab: closet.local:/mnt/nfs /mnt/nfsc/ nfs defaults 0 0
//command = " \!"echo `hostname`>>/mnt/nfsc/hi;cat /mnt/nfsc/hi"; // 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 = Run Program(
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" );
);
pythonCodeFileName = "DeleteMeFrom_" || Regex( Include File List()[1], "([^\!\\!\]*?\.jsl)", "\1" ) || "_.py";
pythonCodeFilePathWin = "\\VBOXSVR\nfsc\farm\";
pythonCodeFilePathLinux = "/mnt/nfsc/farm/";
xdimension = 1920/2;
ydimension = 1080/2;
squaresize = 10; // squaresize^2 processes will run (30 works, slow to launch)
subImageWidth = xdimension / squaresize; // needed to unpack .pbd
if(subImageWidth != floor(subImageWidth),throw("subImageWidth "||char(subImageWidth)));
subImageHeight = ydimension / squaresize; // maybe not needed
if(subImageHeight != floor(subImageHeight),throw("subImageHeight "||char(subImageHeight)));
edges = [=>]; // holds the x/y lef/rig/top/bot mpfr strings for labeling graph xlef
progressbarlength = 50; // in the monitor window
resultWindowScale = 1; // scale preview
digitPrecision = 25; // prec
drillDepth = 800; // deep
colorpower = 1.3;
savetextfile(pythonCodeFilePathWin||pythonCodeFileName,"\[import threading, time, sys # , tty, termios
#import sys, getopt, time
import math
import numpy as np
import gmpy2
from gmpy2 import mpfr
from gmpy2 import mpc
from PIL import Image
eventWorkAvailable = threading.Event()
eventWorkAvailable.clear()
eventWorkCaptured = threading.Event()
eventWorkCaptured.clear()
# all parm names must be exactly 4 chars...simple parser below...
parms = { # all nums kept as strings because the parser is very simple
'xcen': "0.0", # left edge
'xrad': "0.0", # radius, xmax = xmin + 2 * xrad
'ycen': "0.0", # top edge
'yrad': "0.0",
'xpix': "0", # bitmap size
'ypix': "0", # bitmap size
'deep': "0", # depth...max iterations
'xlef': "", # edges of whole pic (not just my bit)
'xrig': "",
'ytop': "",
'ybot': "",
# 'redp': "", # rgb exponents
# 'grnp': "",
# 'blup': "",
# next 3 calculate a new xcen/ycen/xrad/yrad, in place, from current values
'xzfr': "0.0", # for the zoom command, 0..1 of xmin..xmax
'yzfr': "0.0", # for the zoom command, 0..1 of ymin..ymax
'zfra': "0.0", # typ .5 to 2.0, the zoom in/out fraction
'name': "/mnt/xxx/abc001.pbd", # output file name
'ncpu': "100", # number of workers total, must be a square 1,4,16,...100(sweet spot?),...256 may be near a windows limit for puttys
'icpu': "0", # worker id in square grid of ncpu, 0 based, left to right starting at top (give all workers same xcen, etc for whole pic)
'prec': "100", # precision, in base 10 digits (increases with 'deep' and smaller xrad/yrad)
'prog': "." # progressbar? no comma on last one. ***READ ONLY*** progress bar status
}
gmpy2.get_context().precision = int(3.32*int(parms['prec']))
command = "" # global for main worker thread to copy
drop = 0 # global to communicate to main worker thread "drop the picture in progress asap and wait for next command"
def listener(): # listen to stdin for parms and commands
global parms
global command
global drop
global digits
gmpy2.get_context().precision = int(3.32*int(parms['prec'])) # roughly digits # The thread gets its own context
formatstring = '{0:.' + str(int(parms['prec'])+2) + 'G}' # save a few extra, general format
ndrop=0
while 1:
line = sys.stdin.readline().strip()
if not line: # eof detected and turns into normal quit command
line = "quit"
cmd = line[:4] # commands/parms are all exactly 4 chars
op = line[4:].strip() # strip space between cmd...operand and trailing newline
if cmd in parms: # asynchronous set a parm
# cmd should NOT be 'prog', but any other 'xcen'...'prec' parm to be set to a value
parms[cmd] = op # set a parameter value
if cmd == 'prec':
gmpy2.get_context().precision = int(3.32*int(parms['prec'])) # just this thread for now
formatstring = '{0:.' + str(int(parms['prec'])+2) + 'G}' # save a few extra, general format
elif cmd == "read": # asynchronous retrieve a parm, including 'prog' progressbar value
# a leading space...often the first character arrives early and
# gets separated from the message. This gets trimed if present.
print(" "+op,parms[op]) # retrieve a parameter value
sys.stdout.flush()
elif cmd == "drop": # asynchronous
drop = 1 # drop the current picture in progress
ndrop+=1
# parms["prog"]="drp"+str(ndrop)
#
# print(" prog",parms["prog"]) # announce dropped pic
# sys.stdout.flush()
elif cmd == "zoom": # asynchronous
# tiny mpfr calculator for zooming in/out ***using current prec***
xcen = mpfr(parms["xcen"])
xrad = mpfr(parms["xrad"])
xmin = xcen - xrad
xmax = xcen + xrad
ycen = mpfr(parms["ycen"])
yrad = mpfr(parms["yrad"])
ymin = ycen - yrad
ymax = ycen + yrad
xzfr = mpfr(parms["xzfr"])
yzfr = mpfr(parms["yzfr"])
zfra = mpfr(parms["zfra"])
parms["xcen"] = formatstring.format(xmin + xzfr * (xmax - xmin))
parms["ycen"] = formatstring.format(ymin + yzfr * (ymax - ymin))
parms["xrad"] = formatstring.format(xrad * zfra)
parms["yrad"] = formatstring.format(yrad * zfra)
else: # synchronous: these are interpreted below
# quit or make. other garbage is ignored, mostly, but will ***block*** until prev make finishes
command = cmd
#operand = op
eventWorkAvailable.set() # main wakes up and starts capturing details
eventWorkCaptured.wait() # main signals the details are captured (***block*** happens here if make in progress)
eventWorkCaptured.clear() # main will need to set this
if cmd == "quit": # and quit is also interpreted here so it doesn't read again
break
# start the listener thread to process commands on stdin
stdiothread = threading.Thread(target=listener,args=(()))
stdiothread.start()
def mandelbrot(z,maxiter,horizon,log_horizon,log_2):
c = z
for n in range(maxiter):
az = abs(z)
if az > horizon:
return float(n) - float(gmpy2.log(gmpy2.log(az))/log_2 + log_horizon) # https://gist.github.com/jfpuget/16eb22a3c26b6275d9dd -- the smooth log part
z = z*z + c
# return 0
return float(maxiter)
def mandelbrot_set(xmin,xstep,ymin,ystep,nx,ny,maxiter):
global parms
global drop
horizon = mpfr(2.0 ** 40)
log_2 = gmpy2.log(2)
log_horizon = gmpy2.log(gmpy2.log(horizon))/log_2
step10 = min(xstep,ystep)/10 # radius of small triangle of slopes about the point
# power=1
# dx0 = mpfr(step10 * math.sin(0))
# dy0 = mpfr(step10 * math.cos(0))
# dx1 = mpfr(step10 * math.sin(math.pi*120/180))
# dy1 = mpfr(step10 * math.cos(math.pi*120/180))
# dx2 = mpfr(step10 * math.sin(math.pi*240/180))
# dy2 = mpfr(step10 * math.cos(math.pi*240/180))
# old code
# redp = float(parms["redp"])
# grnp = float(parms["grnp"])
# blup = float(parms["blup"])
# n3 = np.empty((ny,nx,3),np.uint8)
# new code
n2 = np.empty((ny,nx),dtype=np.float64)
x = xmin
for ix in range(nx):
#time.sleep(0.001)
if drop == 1:
break
parms["prog"] = str(int(100*ix/nx))+"%"
if ix % 3 == 0:
print(" prog",parms["prog"]) # announce %
sys.stdout.flush()
y = ymin
for jy in range(ny):
if drop == 1:
break
# p0 = mandelbrot( mpc(x+dx0, y+dy0), maxiter, horizon, log_horizon, log_2)
# p1 = mandelbrot( mpc(x+dx1, y+dy1), maxiter, horizon, log_horizon, log_2)
# p2 = mandelbrot( mpc(x+dx2, y+dy2), maxiter, horizon, log_horizon, log_2)
# close = (max(p0,p1,p2)-min(p0,p1,p2))/((p0+p1+p2)/3) # 1,2,3 -> 2/6 501,502,503-> 2/502 0,500,1000 -> 1000/500
# if close > 1:
# close = 1
# close = close
# d = math.atan2(p0,p1)/(math.pi/2) # -1 ... + 1
# e = math.atan2(p1,p2)/(math.pi/2)
# f = math.atan2(p2,p0)/(math.pi/2)
# n3[jy,ix] = (255*close*(d+1)/2,255*close*(e+1)/2,255*close*(f+1)/2)
p = mandelbrot( mpc(x, y), maxiter, horizon, log_horizon, log_2)/maxiter
# old code
# try:
# n3[jy,ix] = (255*p**redp,255*p**grnp,255*p**blup)
# except:
# print(p,type(p),redp,grnp,blup)
# new code
n2[jy,ix] = p
y = y + ystep
x = x + xstep
return n2
while 1: # process work...
eventWorkAvailable.wait()
eventWorkAvailable.clear() # tell the listener we are not available
copycommand = command
command = ""
eventWorkCaptured.set() # tell the listener we've captured the details
if copycommand=="quit":
break
if copycommand=="make":
gmpy2.get_context().precision = int(3.32*int(parms['prec']))
formatstring = '{0:.' + str(int(parms['prec'])+2) + 'G}'
drop = 0 # if this becomes 1, drop this picture and return to waiting
# parms["prog"]="start" # typically 1/8192 step size, this is a smaller, barely started, value
# print(" prog",parms["prog"]) # announce pre-percent
# sys.stdout.flush()
# figure out my responsibility in the grid
mynumber = int(parms["icpu"]) # 0..99
totalnumber = int(parms["ncpu"]) # 100
sqrtnumber = math.sqrt(totalnumber) # 10
myxgridLeft = int(mynumber % sqrtnumber) # 0..9
myygridTop = int(mynumber / sqrtnumber) # 0..9
# make 0..1 fractions of my grid position
myxgridLeft = myxgridLeft / sqrtnumber
myygridTop = myygridTop / sqrtnumber
# unpack the parms for the group's grid responsibility, apply my part
xcen = mpfr(parms["xcen"]) # center and size for the group picture
xrad = mpfr(parms["xrad"])
xmin = xcen - xrad # for the group
parms["xlef"] = formatstring.format(xmin) # record group picture edges
xsize = 2 * xrad # for the group
parms["xrig"] = formatstring.format(xmin+xsize) # record group picture edges
xmin += myxgridLeft * xsize # for me
xsize /= sqrtnumber # for me
ycen = mpfr(parms["ycen"]) # center and size for the group picture
yrad = mpfr(parms["yrad"])
ymin = ycen - yrad # for the group
parms["ytop"] = formatstring.format(ymin) # record group picture edges
ysize = 2 * yrad # for the group
parms["ybot"] = formatstring.format(ymin+ysize) # record group picture edges
ymin += myygridTop * ysize # for me
ysize /= sqrtnumber # for me
# only need to say this once (every thread says it)
print( " xlef",parms["xlef"],
"\nxrig",parms["xrig"],
"\nytop",parms["ytop"],
"\nybot",parms["ybot"],
"\nprec",parms["prec"],
"\ndeep",parms["deep"],
"\nxzfr",parms["xzfr"],
"\nyzfr",parms["yzfr"],
"\nzfra",parms["zfra"],
"\nxcen",parms["xcen"],
"\nxrad",parms["xrad"],
"\nycen",parms["ycen"],
"\nyrad",parms["yrad"],
"\nprog","go go go","\n") # announce edge parms
sys.stdout.flush()
# I'm only responsible for my fraction of the pixels
width = int(int(parms["xpix"]) / sqrtnumber)
height = int(int(parms["ypix"]) / sqrtnumber)
# make the picture
z = mandelbrot_set(xmin,xsize/width,ymin,ysize/height,width,height,int(parms["deep"]))
if drop == 0:
# parms["prog"]="save" # after 100%
# print(" prog",parms["prog"]) # announce saving
# sys.stdout.flush()
# old code returned a png file...
# img = Image.fromarray(z)
# img.save(parms["name"])
# new code now binary... ".pbd" extension "Python Binary Data"...I made that up.
fileobj = open(parms["name"], mode='wb')
z.tofile(fileobj)
fileobj.close()
parms["prog"]="okpic"
else:
parms["prog"]="*BAD*"
print(" prog",parms["prog"]) # announce fail to make pic or good pic
sys.stdout.flush()
# post-quit cleanup
#print("wait for thread to join")
stdiothread.join() # wait until listener actually finishes
print("done bye")
]\");
dtThumb = New Table( "thumbnails",
New Column( "thumb Picture", Expression, "None", Set Display Width( 88 ) ),
New Column( "thumb Time", Format( "hr:m:s", 13, 0 ), Input Format( "hr:m:s", 0 ), Set Display Width( 74 ) ),
New Column( "thumb Data", Expression, "None", Set Display Width( 68 ) ),
New Column( "thumb CenterX", Character, Set Display Width( 82 ) ),
New Column( "thumb RadiusX", Character, Set Display Width( 86 ) ),
New Column( "thumb CenterY", Character, Set Display Width( 78 ) ),
New Column( "thumb RadiusY", Character, Set Display Width( 80 ) ),
New Column( "thumb Left", Character, Set Display Width( 68 ) ),
New Column( "thumb Right", Character, Set Display Width( 67 ) ),
New Column( "thumb Top", Character, Set Display Width( 73 ) ),
New Column( "thumb Bottom", Character, Set Display Width( 84 ) ),
New Column( "thumb Digits", Format( "Fixed Dec", 12, 0 ), Set Display Width( 93 ) ),
New Column( "thumb Iteration", Format( "Fixed Dec", 12, 0 ), Set Display Width( 79 ) )
);
/*
1:aaa
2:bbb
3:dddd
4:eeee
5:fffff
6:gggggg
7:hhhhhh
8:iiiiii
here's a pattern to spatially separate the machines
defaghibfdeghiagfhdeibghfiadeghifbghi
use items from this list, in a ring, to balance the load (aaa is a weaker machine than gggggg)
and distribute the load spatially (don't use aa next to each other in the bitmap, they might be
both in a hard patch)
built by ababab is len 6, put ghi between each is len 24, put de every 6 is len 32, f every 6 is len 37
the prime 37 should work with a 16x16 square OK
*/
separatePlan = [3 4 5 1 6 7 8 2 5 3 4 6 7 8 1 6 5 7 3 4 8 2 6 7 5 8 1 3 4 6 7 8 5 2 6 7 8];
separateIndex=1;
findMachine = Function( {},
{iCondat, n}, // find a machine that is not filled up
separateIndex+=1;
if(separateIndex>nitems(separatePlan),separateIndex=1);
return(separatePlan[separateIndex]);
);
// build/replace/insert new image
histcolor = heatcolor((0::1022)/1022,"spectral")||0;
//shortHigh = j(1,128,1);
//longHigh = j(1,256,1);
//shortLow = j(1,128,0);
//longLow = j(1,256,0);
//shortRiseFromMiddle = (128::255)/255;
//longRiseFromBottom = (0::255)/255;
//shortFallToMiddle = (255::128)/255;
//longFallToBottom = (255::0)/255;
//red = longLow || longLow || longLow || longRiseFromBottom;
//red[1024]=0;
//green = longRiseFromBottom || shortHigh || longHigh || longHigh || shortFallToMiddle;
//green[1024]=0;
//blue = shortRiseFromMiddle || longHigh || longFallToBottom || shortLow || longRiseFromBottom;
//blue[1024]=.25;
//histcolor = rgbcolor( red,green,blue);
generateImage = Function( {},
{minmaxrange, normalized, img, imgseg, hist, color},
minmaxrange = Range( compositionMatrix );
normalized = (compositionMatrix - minmaxrange[1]) / (minmaxrange[2] - minmaxrange[1]);
hist = Sort Ascending( normalized[0] );
hist = hist[1 :: N Rows( hist ) :: N Rows( hist ) / 1024];
loc1 = Loc( hist, 1 );
If( N Items( loc1 ),
hist = hist[1 :: loc1[1]]
); // only one slot for max
posterized = Round( (1024 / N Rows( hist )) * Loc Sorted( hist, normalized ^ colorpower ) ); // <<<<<<<<< trim power
color = Shape( histcolor[posterized], N Rows( compositionMatrix ), N Cols( compositionMatrix ) );
img = New Image( color );
//newwindow("x",img)
// img = newimage("rgb",{normalized^redpower,normalized^greenpower,normalized^bluepower});
Try(
imgseg = resultGB[framebox( 1 )] << FindSeg( PictSeg( 1 ) );
try(imgseg << remove); // delete old and replace with new...fail expected 1st time
resultGB[framebox( 1 )] << Add Image( image( img ), bounds( top( 0 ), Left( 0 ), bottom( 1 ), Right( 1 ) ) );//
imgseg = resultGB[framebox( 1 )] << FindSeg( PictSeg( 1 ) );
imgseg << lock( 1 ); // so the mouse won't interact with the image and drag it around
,//
show(exception_msg)
);
return(img);
);
// launch a thread in the farm
launch = Function( {imagenumber,machine = findMachine()},
{command},
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.8 " || pythonCodeFilePathLinux || pythonCodeFileName;
imgNumToRP[imagenumber] = Run Program(
executable( "plink" ),
options(
localkey || /* -v */ " -hostkey " || condat[machine][hostkey] || " " || remoteuser || condat[machine][server] || " " || command || //
""
),
Parameter( Eval List( {imagenumber, machine} ) ),
readfunction(
Function( {this, parm},
{t, txt, txttrim, imagenumber, machine},
{imagenumber, machine} = parm;
txt = "";
While( this << canread,
t = (this << read);
txt = txt || t;
if(endswith(txt,"\!n")|endswith(txt,"\!r"),break());
wait(.01);
);
/////////////////////
if(contains(txt,"rror"), write(txt));
if(contains(txt,"okpic") & imgNumToIsStoredInGraph[imagenumber] != generation,
// grab the picture from OKPIC
local({left,right,top,bot,filename,rc,imgseg,bindata,matdata,img},
//show(imagenumber);
left = Mod( (imagenumber - 1), squaresize );
col = left * subImageWidth;
right = left + 1;
top = Floor( (imagenumber - 1) / squaresize );
row = top * subImageHeight;
bot = top + 1;
left = left / squaresize;
right = right / squaresize;
top = top / squaresize;
bot = bot / squaresize;
filename=Eval Insert( pythonCodeFilePathWin || "mbrot_^generation^_^Char( imagenumber )^.pbd" );
try(//
//show(filename);
bindata = Load Text File( filename, blob() );
matdata = Blob To Matrix( bindata, "float", 8, "little", subImageWidth );
compositionMatrix[(row+1)::(row+subImageHeight), (col+1)::(col+subImageWidth)] = matdata[0,0];
img = generateImage();// <<<<<<<<<<<<<<< resultGB<<inval;//generateImage();
imgNumToIsStoredInGraph[imagenumber] = generation;
imagesReceived += 1;
if(imagesReceived==squaresize^2,
show(edges);
dtThumb<<addrows(1);
dtThumb:thumbTime[nrows(dtThumb)] = TickSeconds()-StartTime;
dtThumb:thumbPicture[nrows(dtThumb)] = img;
dtThumb:thumbData[nrows(dtThumb)] = compositionMatrix;
dtThumb:thumbIteration[nrows(dtThumb)] = edges["deep"];
dtThumb:thumbDigits[nrows(dtThumb)] = edges["prec"];
dtThumb:thumbCenterX[nrows(dtThumb)] = edges["xcen"];
dtThumb:thumbLeft[nrows(dtThumb)] = edges["xlef"];
dtThumb:thumbRadiusX[nrows(dtThumb)] = edges["xrad"];
dtThumb:thumbRight[nrows(dtThumb)] = edges["xrig"];
dtThumb:thumbBottom[nrows(dtThumb)] = edges["ybot"];
dtThumb:thumbCenterY[nrows(dtThumb)] = edges["ycen"];
dtThumb:thumbRadiusY[nrows(dtThumb)] = edges["yrad"];
dtThumb:thumbTop[nrows(dtThumb)] = edges["ytop"];
)
,//
show(exception_msg);//
);//
);
);
// check all lines for xlef/xrig/ytop/ybot reports and stash them away
gotedge=0;
patmatch(txt,("xlef"|"xrig"|"ytop"|"ybot"|"prec"|"deep"|"xcen"|"xrad"|"ycen"|"yrad")>>edgekey + patspan(" ") +
patbreak("\!n")>>edgevalue + pattest(edges[edgekey]=edgevalue;gotedge=1;1) + "\!n" + patfail(/*get all of them at once*/));
if(gotedge,
titleabove<<settext(edges["ytop"]);
titlebelow<<settext(edges["ybot"]);
titleleft<<settext(edges["xlef"]);
titleright<<settext(edges["xrig"]);
);
txttrim = regex(txt,"([^\!n\!r]+)[\!n\!r]+$","\1"); // keep just the most recent (last) line
//show("trim",chartoblob(txt),ismissing(txttrim),txttrim);
if(!ismissing(txttrim) & length(txttrim) > 5 & try(!contains(imgNumToMonTextBox[imagenumber]<<gettext,"bye"),1),//
Try(//
if(contains(txttrim,"ERROR"),
if(try(imgNumToMachine[imagenumber]<<isreadeof,1),
imgNumToMachine[imagenumber] = 0;
show(txttrim,imagenumber,machine,"relaunching");
launch(imagenumber, machine);
);
txttrim="aaaaa"||char(machine);
);
txttrim = trim(Substr( txttrim, 6 ));//
if(contains(txttrim,"%"),
try(imgNumToMonTextIfBox[imagenumber]<<set(0));
siz=round(progressbarlength*num(txttrim));
try((imgNumToMonProgressBox[imagenumber]<<child)<<width(siz));
try(((imgNumToMonProgressBox[imagenumber]<<child)<<sib)<<width(progressbarlength-siz));
imgNumToMonProgressIfBox[imagenumber]<<set(1);
,
if(txttrim != "." & txttrim != "go go go",//speedup?
//show(txttrim);
try(imgNumToMonTextBox[imagenumber] << settext( left(txttrim,10) ));//
);
try(imgNumToMonTextIfBox[imagenumber]<<set(1));
try(imgNumToMonProgressIfBox[imagenumber]<<set(0));
);
imgNumToResponse[imagenumber] = txttrim;
// imgNumToMonTextBox[imagenumber] << updatewindow;//
, //
Show( exception_msg,"bad", chartoblob(txt) );//
);//
);//
////////////////////////////////////
);
);
);
//Show( imagenumber );
imgNumToRP[imagenumber] << Write( "read prog\!n" );
return( machinename[machine] );
);
imgNumToRP = [=> ]; // holder for active RunProgram objects
imgNumToRPdone = [=> ]; // holder for RunProgram objects that need a moment to terminate
imgNumToMachine = [=> 0]; // which machine is this running on?
imgNumToResponse = [=> ""]; // f(x) labels for graphs
imgNumToMonTextBox = [=> ];
imgNumToMonTextIfBox = [=> ];
imgNumToMonProgressBox = [=> ];
imgNumToMonProgressIfBox = [=> ];
imgNumToIsStoredInGraph = [=>-1];
compositionMatrix = j(ydimension,xdimension,1);//imgNumToSegInGraph = [=>];
stopImgNum = Function( {imgNum},
imgNumToRP[imgNum] << Write( "quit\!n" )
);
stopall = Function( {},
{k, rp, i},
dropall();
For Each( {{k, rp}, i}, imgNumToRP, stopImgNum( k ) );
print("all stopping messages sent");
);
dropall = Function( {},
{k, rp, i},
For Each( {{k, rp}, i}, imgNumToRP, rp << Write( "drop\!n" ) );
);
startdrag=function({x,y},
if(ismissing(dragpos0[1]),
dragpos0[1]=x;dragpos0[2]=y;//
dragpos1[1]=x;dragpos1[2]=y;//
// show("cap",dragpos0);
,//
dragpos1[1]=x;dragpos1[2]=y;//
// show("-->",dragpos1);
)
);
enddrag=function({x,y},
dragpos1[1]=x;dragpos1[2]=y;//
dropall(); // stop old work
/*
redp ^redpower^
grnp ^greenpower^
blup ^bluepower^
*/
commands =
"prec ^digitPrecision^
deep ^drillDepth^
xzfr ^(dragpos0[1]+dragpos1[1])/2^
yzfr ^(dragpos0[2]+dragpos1[2])/2^
zfra ^max(abs(dragpos0[1]-dragpos1[1]),abs(dragpos0[2]-dragpos1[2]))^
zoom
name ^pythonCodeFilePathLinux^mbrot_^generation^_^Char( k )^.pbd
make
";
imgNumToIsStoredInGraph = [=>-1];
//foreach({{k,v},i},imgNumToSegInGraph,v<<remove;);
compositionMatrix[0,0] = 1; //imgNumToSegInGraph = [=>];
imagesReceived = 0;
sendCommandsToAll_k_rp_i(commands);
dragpos0[1]=.;dragpos0[2]=.;
dragpos1[1]=.;dragpos1[2]=.;
StartTime=tickseconds();
);
imagesReceived = 0;
StartTime=tickseconds();
dragpos0=[. .];
dragpos1=[.,.];
resultWindow = New Window( "result",
V List Box(
H Center Box( titleabove = Text Box( "above" ) ),
H List Box(
V Center Box( titleleft = Text Box( "left", <<Rotate Text( "left" ) ) ),
resultGB = Graph Box(//
framesize( xdimension * resultWindowScale, ydimension * resultWindowScale ),
X Scale( 0, 1 ),
suppressaxes,
<<yaxis( Min( 1 ), Max( 0 ) ),
<<AddGraphicsScript("front",
mousetrap(startdrag(x,y),enddrag(x,y));
//nope. not here. generateImage(); // can this be here? it might be modifying the display tree during show...
if( !ismissing(dragpos0[1]),
pencolor("red");
linestyle("dotted");
line({dragpos0[1],dragpos0[2]},{dragpos0[1],y});
line({dragpos0[1],y},{x,y});
line({x,y},{x,dragpos0[2]});
line({x,dragpos0[2]},{dragpos0[1],dragpos0[2]});
//rect(min(dragpos0[1],x),max(dragpos0[2],y),max(dragpos0[1],x),min(dragpos0[2],y),1)
);
)
),
V Center Box( titleright = Text Box( "right", <<Rotate Text( "left" ) ) )
),
H Center Box( titlebelow = Text Box( "below" ) ),
hlistbox(
textbox("color"),sliderbox(.1,10,colorpower,generateimage()),
spacerbox(size(20,1)),
globalbox(drillDepth),
spacerbox(size(20,1)),
globalbox(digitPrecision)
)
)
);
monitorWindow = New Window( "monitor",
V List Box(//
H List Box(//
Button Box( "stop all", stopall() ), //
// Button Box( "drop all", dropall() ) //
), //
lub = Lineup Box( N Col( squaresize ),
spacing( 3, 3 ),
<<backgroundcolor(rgbcolor(.94,.97,.97))
) //
) //
);
write("\!ncheckpoint 1\!n");wait(.1);
monbox = Function( {imgNum, macnam},
Eval(
Eval Expr(
H List Box(
Spacer Box( size( 1, 30 ), color( rgbcolor(1,1,1) ) ),// left edge
Spacer Box( size( 5, 1 ), color( rgbcolor(1,1,1) ) ),// first part of top edge
V List Box(
Spacer Box( size( progressbarlength, 1 ), color( rgbcolor(1,1,1) ) ), // most of top edge
spacerbox(size(1,5)),
imgNumToMonTextIfBox[imgNum] = If Box( 1, imgNumToMonTextBox[imgNum] = Text Box( Char( imgNum ) ) ),
imgNumToMonProgressIfBox[imgNum] = If Box(
0,
imgNumToMonProgressBox[imgNum] = H List Box(
Spacer Box( size( 0, 7 ), color( "green" ) ),
Spacer Box( size( progressbarlength, 7 ), color( rgbcolor(.9,.9,.9) ) )
)
),
Text Box( macnam )
),
<<backgroundcolor(rgbcolor(.95,.95,.95))
)
)
)
);
write("\!ncheckpoint 2\!n");wait(.1);
wait(.1);monitorWindow<<movewindow(0,0);wait(.1);
write("\!ncheckpoint 3\!n");wait(.1);
For( i = 1, i <= squaresize^2, i += 1,
lub << append( monbox( i, launch( i ) ) );
//wait(0);
);
write("\!ncheckpoint 4\!n");wait(.1);
generation = 0;
sendCommandsToAll_k_rp_i = function({commands},
if(contains(commands,"make"), generation+=1;show(generation));
imgNumToResponse = [=> ""];
For Each( {{k, rp},i}, imgNumToRP,
rp << Write( regex(evalinsert(commands),"[\!r\!n]+","\!n",GLOBALREPLACE ));// change returns to newlines for linux
);
show("broadcast commands done");
);
// bring it up. the threads need some coaxing.
commands =
"zfra ^z^
read zfra
";
write("\!ncheckpoint 5\!n");wait(.1);
For( z = 1, z < 10, z += 1,
sendCommandsToAll_k_rp_i( Eval Insert( commands ) );
Print( "Ready", z );
start = tickseconds();
While( N Items( imgNumToResponse ) < N Items( imgNumToRP ),
if(tickseconds()-start > 10 & N Items( imgNumToResponse ) == 0, throw("nothing is happening"));
write(evalinsert("\!nwaiting for responces ^N Items( imgNumToResponse )^, ^N Items( imgNumToRP )^"));
Wait( 1 );
);
needmore = 0;
For Each( {{k, v}}, imgNumToResponse, If( trim(v) != char(z), needmore = 1 ) );
if(needmore==0,break());
);
write("\!ncheckpoint 6\!n");wait(.1);
begin=tickseconds();
// make same image... a 10x10 processor grid on the beowulf cluster seems good.
// 1080 /2/2/2/3/3/3/5
// 1920 /2/2/2/2/2/2/2/3/5
// 2 3 4 5 6 8 10 12 15 20 24 30 40
/*
x -1.94062988168240384898793330252067703470368414336
x -1.94062988168240384898793330252067703470368414330
y 0.00100077622840272515587559395423166491124004188
y 0.00100077622840272515587559395423166491124004184
*/
/*
redp ^redpower^
grnp ^greenpower^
blup ^bluepower^
xcen -0.5396920368085179367059520359131710652013065
ycen -0.6611804024210354957147719385629364551294705
ycen -0.049842152966738802163856592043818211312978962993082255019266085623
-0.049842152966738802163856592043818211312978962993082255019266085622
xcen 0.336447185016324548337106220487657679616547644208117764858297769859
0.336447185016324548337106220487657679616547644208117764858297769858
*/
commands =
"read prog
prec ^digitPrecision^
deep ^drillDepth^
xcen -1.9406298816824038489879333025206770347036841433
xrad 3
ycen 0.0010007762284027251558755939542316649112400418
yrad 2
xpix ^xdimension^
ypix ^ydimension^
ncpu ^squaresize*squaresize^
icpu ^i-1^
name ^pythonCodeFilePathLinux^mbrot_^generation^_^Char( k )^.pbd
make
";
sendCommandsToAll_k_rp_i(commands);
print("polling...");// check for okpic and *BAD*
keepGoing = 1;
While( keepGoing,
keepGoing = 0;
For Each( {{k, rp},i}, imgNumToRP,
if( (rp<<CanRead) | (rp<<CanWrite()) | !(rp<<IsReadEOF),
if( trim(imgNumToResponse[k]) != "okpic" & trim(imgNumToResponse[k]) != "*BAD*",
keepGoing = 1;//
commands = "read prog\!n";
rp << Write( commands );
wait(1);//
);
);
);
);
show(tickseconds()-begin,edges);
print("JSL stopped");
- © 2024 JMP Statistical Discovery LLC. All Rights Reserved.
- Terms of Use
- Privacy Statement
- About JMP
- JMP Software
- JMP User Community
- Contact
You must be a registered user to add a comment. If you've already registered, sign in. Otherwise, register and sign in.