Sunday, 16 November 2014

A man is rich in proportion to the number of things he can afford to let alone.

"When he has succeeded in being Caesar's friend, he has none the less failed to get what he was seeking. For what is it that every man is seeking? To live securely, to be happy, to do everything as he wishes to do, not to be hindered, not to be subject to compulsion … If he is not invited [by Caesar], he is hurt, and if he is invited, he dines like a slave at a master's table, all the time careful not to say or do something foolish … When did you sleep more peacefully, now or before you became Caesar's friend?" Epictetus. Book IV, Chapter 1.

"Men in great place are thrice servants: servants to the sovereign or state, servants of fame, and servants of business, so as they have no freedom, neither in their persons, nor in their actions, nor in their times." Francis Bacon. Essays, Of Great Place.

"When, therefore, neither those who are styled kings live as they will, nor the friends of these kings, what free men are left? Seek and you will find. For nature has given you resources to find the truth." Epictetus. Book IV, Chapter 1.

Wednesday, 8 October 2014

Trippy plots and codes...

I noticed wikipedia has a new visualization for the Riemann zeta function. However I'm not sure it's quite best way to show it in every case. This is how it looks:

This is a complex function so each point in 2D-plane has phase and amplitude. In the domain coloring scheme the phase is represented by hue (color) and amplitude by value (intensity). This particular image has some modulo coloring of the amplitude etc. and the difference can be most easily seen by comparing it to these 3D-plots I made where the hue independently represents the phase and amplitude is simply represented by the 3rd dimension.

The second image has phase exaggerated quite a lot to make it clear that for these values the true change of phase is quite different from what it appears in the wikipedia image.

While I was at it, I also plotted the Mandelbrot fractal in terms of the complex value after iteration instead of number of iterations as it's usually done. This too resulted in a rather trippy image. Though obviously hue coloring tends to always results in rather trippy images.

Here hue represents phase angle of the complex number and value (intensity) the amplitude. 3D representation is also rather interesting.
Here only hue is varied.

N = 1024;
s = zeros(N, N);
for y = 1:N
    for x = 1:N
        z0 = (6*x/N-3.5) + j*(6*y/N-3);
        z = 0;
        iteration = 0;
        max_iteration = 256;
        while norm(z) < 8 && iteration < max_iteration
            z = z^2 + z0;
            iteration = iteration + 1;
        s(x, y) = log(iteration);
        S(x, y) = z;


h = angle(S);
h = h-min(min(h));
h = h/max(max(h));

v = abs(S);

hsv(:,:,1) = h;
hsv(:,:,2) = 1;
hsv(:,:,3) = 1;
% hsv(:,:,3) = v/max(max(v));
rgb = hsv2rgb(hsv);

x = 6*[1:N]/N-3.5;
y = 6*[1:N]/N-3;
surf(x, y, v, rgb, 'edgecolor', 'none');
shading interp;

image(x, y, rgb);
imwrite(rgb, 'compmand.png', 'png')

How to test shaders with minimal effort on OSX (though it'll work with relatively minor changes on linux as well)...

import ctypes
from ctypes import *

d = cdll.LoadLibrary("SceneKit.framework/SceneKit")

libSDL = cdll.LoadLibrary("libSDL.dylib")

a = c_char_p("""
uniform float t;
vec4 orb;
float ss;
float map(vec3 p){
 float scale = 1.0;
 orb = vec4(1000.0);
 for( int i=0; i<8;i++ ){
  p = -1.0 + 2.0*fract(0.5*p+0.5);
  float r2 = dot(p,p);
         orb = min( orb, vec4(abs(p),r2) );
  float k = max(ss/r2,0.1);
  p *= k;
  scale *= k;
 return 0.25*abs(p.y)/scale;
float trace( in vec3 ro, in vec3 rd ){
 float maxd = 100.0;
 float precis = 0.001;
 float h=precis*2.0;
   float t = 0.0;
 for(int i=0; i<200; i++){
  if( abs(h)<precis||t>maxd ) continue;
    t += h;
      h = map( ro+rd*t );
     if( t>maxd ) t=-1.0;
 return t;
vec3 calcNormal( in vec3 pos ){
 vec3  eps = vec3(.0001,0.0,0.0);
 vec3 nor;
 nor.x = map(pos+eps.xyy) - map(pos-eps.xyy);
 nor.y = map(pos+eps.yxy) - map(pos-eps.yxy);
 nor.z = map(pos+eps.yyx) - map(pos-eps.yyx);
 return normalize(nor);
void main(void){
 vec2 iResolution=vec2(512,384);
 vec2 p = -1.0 + 2.0*gl_FragCoord.xy / iResolution.xy;
    p.x *= iResolution.x/iResolution.y;
 float time = t*0.25;
    ss = 1.1 + 0.5*smoothstep( -0.3, 0.3, cos(0.1*t) );
 vec3 ro = vec3( 2.8*cos(0.1+.33*time), 0.4 + 0.30*cos(0.37*time),2.8*cos(0.5+0.35*time) );
 vec3 ta = vec3( 1.9*cos(1.2+.41*time), 0.4 + 0.10*cos(0.27*time), 1.9*cos(2.0+0.38*time) );
 float roll = 0.2*cos(0.1*time);
 vec3 cw = normalize(ta-ro);
 vec3 cp = vec3(sin(roll), cos(roll),0.0);
 vec3 cu = normalize(cross(cw,cp));
 vec3 cv = normalize(cross(cu,cw));
 vec3 rd = normalize( p.x*cu + p.y*cv + 2.0*cw );
 vec3 col = vec3(0.0);
 float t = trace( ro, rd );
 if( t>0.0 ){
  vec4 tra = orb;
  vec3 pos = ro + t*rd;
  vec3 nor = calcNormal( pos );
         vec3  light1 = vec3(  0.577, 0.577, -0.577 );
         vec3  light2 = vec3( -0.707, 0.000,  0.707 );
  float key = clamp( dot( light1, nor ), 0.0, 1.0 );
  float bac = clamp( 0.2 + 0.8*dot( light2, nor ), 0.0, 1.0 );
  float amb = (0.7+0.3*nor.y);
  float ao = pow( clamp(tra.w*2.0,0.0,1.0), 1.2 );
  vec3 brdf  = 1.0*vec3(0.40,0.40,0.40)*amb*ao;
  brdf += 1.0*vec3(1.00,1.00,1.00)*key*ao;
  brdf += 1.0*vec3(0.40,0.40,0.40)*bac*ao;
  vec3 rgb = vec3(1.0);
  rgb = mix( rgb, vec3(1.0,0.80,0.2), clamp(6.0*tra.y,0.0,1.0) );
  rgb = mix( rgb, vec3(1.0,0.55,0.0), pow(clamp(1.0-2.0*tra.z,0.0,1.0),8.0) );
  col = rgb*brdf*exp(-0.2*t);
 col = sqrt(col);
 col = mix( col, smoothstep( 0.0, 1.0, col ), 0.25 );

b = d.glCreateShader(35632)
c = d.glCreateProgram();
d.glShaderSource(b, 1, c_long(ctypes.addressof(a)), 0)
d.glAttachShader(c, b);
b = d.glGetUniformLocation(c, "t");

while 1:
 d.glUniform1f(b, c_float(t));
 d.glRecti(-1, -1, 1, 1);

How about this shit (
[bits 32]
; Use Python to load the SDL library, open an SDL window, mmap us into memory,
; then start the intro.

db `python -c"from ctypes import *;S=CDLL('');`
db `CFUNCTYPE(int)(S.mmap(0,1,7,2,'$0',2),0)+167)`
db `(S.SDL_SetVideoMode(640,480,0,0),S.SDL_PollEvent,S.SDL_Flip);" #`
; Global constants
%define XRES          640        ; X resolution of screen
%define YRES          480        ; Y resolution of screen
%define VOXEL_SIZE    4          ; Size of a voxel (2^VOXEL_SIZE)
%define MAXITS        32         ; Maximum number of iterations when marching
%define ROT_SPEED     120        ; Speed of rotation (larger = slower)

; Pseudonyms for stack addresses within main()

%define SURFACE       ebp + 08 ; Pointer to an SDL_Surface
%define SDL_POLLEVENT ebp + 12 ; Pointer to SDL_Pollevent() function
%define SDL_FLIP      ebp + 16 ; Pointer to SDL_Flip() function

    push ebp
    mov  ebp, esp

This fellow does interesting trick to avoid having to make a proper binary on linux, resulting in 512 byte demo.

One can dynamically (or sometimes even not so dynamically) load function addresses in C as well.

#include <dlfcn.h>
typedef void func(int,int,int,int);
        dlopen("SceneKit.framework/SceneKit", RTLD_LAZY);

        void *a = dlsym(RTLD_DEFAULT, "glCreateShader");
        void *b = dlsym(RTLD_DEFAULT, "glCreateProgram");
        void *c = dlsym(RTLD_DEFAULT, "glShaderSource");
        void *d = dlsym(RTLD_DEFAULT, "glCompileShader");
        void *e = dlsym(RTLD_DEFAULT, "glAttachShader");
        void *f = dlsym(RTLD_DEFAULT, "glLinkProgram");
        void *g = dlsym(RTLD_DEFAULT, "glUseProgram");
        void *h = dlsym(RTLD_DEFAULT, "glGetUniformLocation");
        void *i = dlsym(RTLD_DEFAULT, "glUniform1fv");
        void *j = dlsym(RTLD_DEFAULT, "glRecti");
        void *k = dlsym(RTLD_DEFAULT, "printf");

        printf("glCreateShader: %p\n", a);
        printf("glCreateProgram: %p\n", b);
        printf("glShaderSource: %p\n", c);
        printf("glCompileShader: %p\n", d);
        printf("glAttachShader: %p\n", e);
        printf("glLinkProgram: %p\n", f);
        printf("glUseProgram: %p\n", g);
        printf("glGetUniformLocation: %p\n", h);
        printf("glUniform1fv: %p\n", i);
        printf("glRecti: %p\n", j);
        printf("printf: %p\n", k);

 func* fu = (func*)0x10f7f8619;
 printf("%p\n", &printf);
// f(640,480,0,0);
One can also mmap raw binary from a file quite easily and execute it:
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <sys/mman.h>

 char buffer[1000];

 FILE *f;
 f = fopen("h64", "rb");
 fseek(f, 0, SEEK_END);
 int size = ftell(f);
 fseek(f, 0, SEEK_SET);
 fread(buffer, size, 1, f);

 printf("Code Size: %d\n", size);
 printf("Prot: %d\n", PROT_EXEC|PROT_READ|PROT_WRITE);
 printf("Flags: %d\n", MAP_PRIVATE|MAP_ANON);

 void* mem_map = mmap(0, sizeof(buffer), PROT_EXEC|PROT_READ|PROT_WRITE, MAP_PRIVATE|MAP_ANON, -1, 0);
 printf("Address: %p\n", mem_map);
 memcpy(mem_map, buffer, sizeof(buffer));


...and the asm (this is for MacOS so if you want to do it for linux, alter it a bit):

; nasm -f bin binary64.asm
[BITS 64]

 mov rax, 0x2000004  ; write
 mov rdi, 1   ; stdout
 lea rsi, [rel msg]
 mov rdx, msg.len


msg: db 'Hello!',0x0A
.len: equ $ - msg

Native OpenGL on Mac OS X:
python -c"from ctypes import *;from Cocoa import *;a=c_char_p(\"float t=gl_Color.x*999.;int m=6;float c=1./(3.*pow(2.,float(m)));vec2 rot(vec2 u,float a){return vec2(u.x*cos(a)-u.y*sin(a),u.y*cos(a)+u.x*sin(a));}void main(void){vec2 u=vec2(640,480);u=-.5*(u-2.*gl_FragCoord.xy)/u.x;u=rot(u,t);u*=sin(t)*.5+1.5;float s=.3;for(int i=0;i<m;i++){u=abs(u)-s;u=rot(u,t);s=s/2.1;}float c=length(u)>c?0.:1.;gl_FragColor=vec4(c,c,c,1.);}\");b=cdll.LoadLibrary(\"SceneKit.framework/SceneKit\");NSApplicationLoad();w=NSWindow.alloc().initWithContentRect_styleMask_backing_defer_(NSMakeRect(0,0,640,480),1,2,1);w.setContentView_(NSOpenGLView.alloc().initWithFrame_pixelFormat_(w.contentRectForFrameRect_(w.frame()),NSOpenGLPixelFormat.alloc().initWithAttributes_([73,5])));w.orderFrontRegardless();d=b.glCreateShader(35632);c=b.glCreateProgram();b.glShaderSource(d,1,c_long(addressof(a)),0);b.glCompileShader(d);b.glAttachShader(c,d);b.glLinkProgram(c);b.glUseProgram(c);
while 1:c+=1;b.glColor3us(c_int(c),0,0);b.glRecti(-1,-1,1,1);b.glSwapAPPLE()"
Looks like one can save 4 bytes in the gzip header technique by just leaving last 4 out. They are crc-32, so they aren't really needed. It'll say truncated input, but works just fine. Previous code gives 675 byte executable (Mac).

echo 'cp \x240 a;sed 1d \x240|zcat>a;./a' > intro
gzip -cn9 > intro1.tmp
dd if=./intro1.tmp of=./intro0.tmp bs=$$((`ls -l intro1.tmp|awk '{ print $$5 }'`-4)) count=1
cat intro0.tmp >> intro
chmod +x intro
rm *.tmp
...\x24 is used in makefile where $ is reserved.

MMAP hello, world in python (binary with nasm):

[bits 32]

db `python -c"from ctypes import *;S=CDLL('');`
db `CFUNCTYPE(int)(S.mmap(0,1,7,2,'$0',2),0)+109)();" #`

 jmp     msg

        mov     eax, 4
        mov     ebx, 1
        pop     ecx
        mov     edx, 7
        int     0x80

msg:    call    start
        db 'Hello!',0xA

Execute binary data in C:

char* code = "\xeb\x13\xb8\x04\x00\x00\x00\xbb\x01\x00\x00\x00\x59\xba\x07\x00\x00\x00\xcd\x80\xc3\xe8\xe8\xff\xff\xff\x48\x65\x6c\x6c\x6f\x21\x0a";

main() {
...stupid blogspot fails to format the text again, but I can't be bothered to give a shit right now.

Monday, 28 April 2014

What if it is all just boring?

I'm not troubled by the fact that there are (mathematical) problems which I do not understand. However, I'm deeply troubled that there are problems which I do understand and yet remain unable to solve.

Suppose the Universe is simply a big computer, a deterministic binary calculator, simple naive information, bit flips and other trivial operations. A machine which simply must exist like the value of pi must be whatever it is, something as simple (in relative terms) as natural numbers for some reason. Suppose we’ll never exceed the speed of light or build a time machine because the Universe doesn’t support those ideas. Suppose all known mysteries of math will have a straightforward closure. P=NP, co-primes can be factored in polynomial time, a perfect games of chess can be found, go gets solved. What goes up, must come down. All fundamental principles will be known. Standard model, quantum mechanics and relativity will merge to give us the theory of everything. Life is just a matter of complexity, mind nothing more than a complex machine, free will an illusion. Life after death doesn’t exit. Art is a simple pattern corresponding to our similarity, brain states, math, and evolutionary history of our species. Love just chemistry. Moore’s law hits a limit. Peace on earth, a cure for cancer. All predictable. All controllable by a perfect machine. AI will surpass human capacity in every way imaginable, but will reveal nothing fundamentally unexpected or unpredictable. Suppose normal biological death will be cured so we may live forever if we so desire and choose.

Then what?

What reason will there be to live? Other than perhaps the human compulsion, the need to feel, but will it be any different from the chemical high of heroin? Was it ever different? Ever noble? Does it need to be? Was it always arbitrary? When all is said and done, would it not be better to remain small, to remain stupid and weak, out of control, be able to forget so we might feel again like the first time? But is this not self deception? Does it have value? Does anything have value? Perhaps not really, but it may not matter, because we are small today and will undoubtedly remain small for years to come. Only in the far distant future can we hope to see far enough on our path to give meaningful answers to these questions.

What about living (forever)? Why does it have value? Why does continuity have value? Because I feel so? But why do I feel so? What difference does it make? Perhaps nothing, and it might be that living a million years would be enough, I wouldn’t know, I can't know, until I’ve lived a million years.

Perhaps you really can’t get an ought from an is. Perhaps it all is what it is, trivial.

I still wish I could have my millions of euros, millions of years and millions of lives, for some reason. I wish to be, do and will. And who’s to say there’s anything wrong with that.

But not to give up on hope entirely... There are still wondrous mysteries and things like the infinity appearing in all of math and science, and many wonderful problems, something which might imply the Universe is still a bit more complicated and wonderful than we can imagine. The quantum mechanics possibly implying there is a reality hidden from us containing an amount of information which is possibly infinitely greater than that which we can observe now (or not). Perhaps there is hope and the Universe is still much greater than we've ever been able to imagine. I just hope I get to experience someday at least some of that greatness.

Halting problem? What problem? It is well know that “...any finite-state machine, if left completely to itself, will fall eventually into a perfectly periodic repetitive pattern. The duration of this repeating pattern cannot exceed the number of internal states of the machine…”

Thursday, 27 February 2014

If I say that I do what I want, it doesn't mean that everything I do is somehow unconventional

An ant starts to crawl along a taut rubber rope 1 km long at a speed of 1 cm per second (relative to the rubber it is crawling on). At the same time, the rope starts to stretch by 1 km per second (both in front of and behind the ant, so that after 1 second it is 2 km long, after 2 seconds it is 3 km long, etc).

When will the ant reach the end of the rope?

...the answer in the end.


Intel Ivy Bridge CPUs (rdrand instruction) and Raspberry Pi have a hardware random number generator which allows the production of seemingly very secure non-deterministic random numbers from thermal etc. sources for cryptographic purposes.

I also noticed that my MacBook Air uses Ivy Bridge i5 which comes with the rdrand instrunction. However the gcc on MacOS didn't seem to support this instruction so I wrote an assembly program (using nasm) to test if I get random numbers out of it. Seemed to work.

I tested the hwrng on the rpi and found that FIPS 140-2 gives similar results to /dev/urandom on linux, however /dev/hwrng on rpi should be intrincically non-deterministic, or at least similarly non-deterministic as thermal noise. The performance seems to be around 107 kB/s, rdrand appeared to give around 240 MB/s on a single core and single thread. I didn't bother tweaking. The official figure is supposedly at least 500 MB/s.

It's another question whether one can trust these hardware implementations or not as we have no sure way of knowing what they actually do and hardware manufacturer may be collaborating with NSA for example and get pseudorandom numbers anyway which might have a pattern identifiable by someone who knows how it was done.

; nasm -f macho64 rdr.asm
; nasm -f elf64 rdr.asm
bits 64
global rdr
section .text
        rdrand rax
;        jnc rdr ; if CF=0, random data not available, try again, however, does not occur on Ivy Bridge so can be ignored at the moment

; gcc rdr.c rdr.o
#include <stdio.h>

int rdr();
int main() {
        int x;
        for(x=0; x<256*256*256*2; x++)
                printf("%c", rdr());

dd if=/dev/urandom of=./urandom.bin count=65536
./a.out > hwrandom.bin
echo --rdrand--rdrand--rdrand--rdrand--
cat hwrandom.bin | rngtest
echo --urandom--urandom--urandom--urandom--
cat urandom.bin | rngtest

rngtest: starting FIPS tests...
rngtest: entropy source exhausted!
rngtest: bits received from input: 268435456
rngtest: FIPS 140-2 successes: 13413
rngtest: FIPS 140-2 failures: 8
rngtest: FIPS 140-2(2001-10-10) Monobit: 1
rngtest: FIPS 140-2(2001-10-10) Poker: 0
rngtest: FIPS 140-2(2001-10-10) Runs: 5
rngtest: FIPS 140-2(2001-10-10) Long run: 2
rngtest: FIPS 140-2(2001-10-10) Continuous run: 0

rdrand: FIPS 140-2 failures: 14, 10, 14, 11, 9, 3, 7, 18, 13, 10
urandom: FIPS 140-2 failures: 11, 14, 9, 9, 7, 10, 10, 9, 13, 11

One should notice that it's supposed to fail occasionally on a perfect random number generator so nothing unusual here.

I ran the tests on one of my servers instead.

$ cat /proc/cpuinfo | grep -i 'model name' | head -n 1
model name      : Intel(R) Core(TM) i5-3570K CPU @ 3.40GHz

$ cat /proc/cpuinfo | grep -i rdrand | head -n 1
flags           : fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush dts acpi mmx fxsr sse sse2 ss ht tm pbe syscall nx rdtscp lm constant_tsc arch_perfmon pebs bts rep_good nopl xtopology nonstop_tsc aperfmp
erf pni pclmulqdq dtes64 monitor ds_cpl vmx est tm2 ssse3 cx16 xtpr pdcm pcid sse4_1 sse4_2 popcnt tsc_deadline_timer aes xsave avx f16c rdrand lahf_lm ida arat epb xsaveopt pln pts dtherm tpr_shadow vnmi flexpriority ept vpid fsgsbase smep erms

It gives "Illegal instruction (core dumped)" if you try to run it on something other than Ivy Bridge.

Similar tests on the rpi...

$ sudo cat /dev/hwrng | rngtest -c 1024
rngtest 2-unofficial-mt.14
Copyright (c) 2004 by Henrique de Moraes Holschuh
This is free software; see the source for copying conditions.  There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

rngtest: starting FIPS tests...
rngtest: bits received from input: 20480032
rngtest: FIPS 140-2 successes: 1023
rngtest: FIPS 140-2 failures: 1
rngtest: FIPS 140-2(2001-10-10) Monobit: 1
rngtest: FIPS 140-2(2001-10-10) Poker: 0
rngtest: FIPS 140-2(2001-10-10) Runs: 0
rngtest: FIPS 140-2(2001-10-10) Long run: 0
rngtest: FIPS 140-2(2001-10-10) Continuous run: 0
rngtest: input channel speed: (min=70.080; avg=958.421; max=1627604.167)Kibits/s
rngtest: FIPS tests speed: (min=841.647; avg=3216.511; max=6401.590)Kibits/s
rngtest: Program run time: 27469963 microseconds

$ sudo cat /dev/urandom | rngtest -c 1024
rngtest 2-unofficial-mt.14
Copyright (c) 2004 by Henrique de Moraes Holschuh
This is free software; see the source for copying conditions.  There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

rngtest: starting FIPS tests...
rngtest: bits received from input: 20480032
rngtest: FIPS 140-2 successes: 1023
rngtest: FIPS 140-2 failures: 1
rngtest: FIPS 140-2(2001-10-10) Monobit: 0
rngtest: FIPS 140-2(2001-10-10) Poker: 0
rngtest: FIPS 140-2(2001-10-10) Runs: 1
rngtest: FIPS 140-2(2001-10-10) Long run: 0
rngtest: FIPS 140-2(2001-10-10) Continuous run: 0
rngtest: input channel speed: (min=1.510; avg=29.660; max=1733.953)Mibits/s
rngtest: FIPS tests speed: (min=843.719; avg=4117.551; max=6403.689)Kibits/s
rngtest: Program run time: 5620869 microseconds

I didn't run the stuff on /dev/random because it's a very low throughput entropy source.

Easiest way to fill the entropy pool (hardware noise):

$ cat /proc/sys/kernel/random/entropy_avail
$ ls -R /
$ cat /proc/sys/kernel/random/entropy_avail

-- -- -- --

GPU accelerated fft on rpi

sudo ./hello_fft.bin  17 8
rel_rms_err = 7e-05, usecs = 17168, k = 0

... this is about 60MS/s, for the xcorr calculation discussed before it would be 2x20 MS/s (2 forwards, 1 backwards).


Eating one banana gives you higher dose of radiation than living next to a nuclear power plant for a year. Sleeping next to another human for 100 days equals one dental x-ray.


…the answer to the ant problem being about 8.8e43421 years.

(A huge, yet finite number much larger than any quantity in this universe.)