Spacetime tells matter how to move, matter tells spacetime how to curve.

"Beyond the corridor of our spacetime there are infinite number of universes, each of them governed by their own set of laws."

Not general relativity, maybe next time if my computer has enough memory and time for that.

That's all Folks!

% ffmpeg -r 30 -f image2 -i %04d.png -ss 00:00:5 -vcodec libx264 -crf 2 output.mp4
clear all; close all; map = parula(256);
N    = 256;
phi  = zeros(N, N);
dphi = zeros(N, N);
ddphi= zeros(N, N);
rho  = zeros(N, N);
[x y] = meshgrid(linspace(-1, 1, N));
f = 0;
for t = 1:3*N
    A = 0.05;
    if t>N
        f = 0.01;
    end
    x1 = cos(2*pi*f*t)*A;
    y1 = sin(2*pi*f*t)*A;
    x2 = -cos(2*pi*f*t)*A;
    y2 = -sin(2*pi*f*t)*A;
    rho = exp(-2000*((x+x1).^2+(y+y1).^2)) + exp(-2000*((x+x2).^2+(y+y2).^2));

    ddphi = rho + del2(phi);
    dphi = dphi + ddphi;
    
     phi(1, :)   =   phi(2, :);
    dphi(1, :)   = -dphi(2, :);
     phi(end, :) =   phi(end-1, :);
    dphi(end, :) = -dphi(end-1, :);
     phi(:, 1)   =   phi(:, 2);
    dphi(:, 1)   = -dphi(:, 2);
     phi(:, end) =   phi(:, end-1);
    dphi(:, end) = -dphi(:, end-1);

    phi = phi + dphi;
    
    cla;
    c = phi/max(max(phi));
    surface(c, 'edgecolor', 'none');
    colormap(map);
    caxis([0 1]);
    axis([1 N 1 N]);
    axis square;
    axis off;
    shading interp;
    drawnow;
    imwrite(c*255, map, ['png/' num2str(t, '%04.f') '.png']);
end

50 Things I'd like to do if I had all the time in the world

 1. Solve quantum gravity.
 2. Crack RSA (factoring) or prove it impossible (polynomial time).
 3. Build a working quantum computer (capable of factoring RSA) or prove it impossible (a separate question from the possibility of classical factoring).
 4. Learn to play piano (I can a bit, but learn it well).
 5. Build my own CPU from transistor level up.
 6. Build my own compiler.
 7. Build my own operating system.
 8. Build my own AI that rivals humans in general creative tasks.
 9. Build my own house.
10. Build my own car.
11. Acquire a few airplanes.
12. Build my own spaceship.
13. Build my own research satellite.
14. Visit the moon.
15. P vs. NP.
16. Climb Aconcagua.
17. Climb Mt. Denali.
18. Climb Mt. Everest.
19. Learn to fly a plane.
20. Visit Mt. Thor (Baffin island).
21. Sail to Ball's Pyramid from Sydney.
22. Visit Svalbard.
23. Directly observe the event horizon of the Milky Way central black hole.
24. Solve the mystery of baryon asymmetry.
25. Solve the mystery of dark energy and dark matter.
26. Crack FTL/time travel or prove them impossible.
27. Find extraterrestrial life.
28. Learn creative drawing.
29. Walk across the Canadian wilderness (~5000 km).
30. If there is a mystery related to human consciousness, solve it, if not, make it obvious for everone there is none.
31. See the day when religions basically are no more.
32. See the day when aging and all disease are cured.
33. See the Milky Way-Andromeda collision.
34. Visit a planet inhabited by estraterrestrial life.
35. Develop more meaningful social relationships.
36. Acquire a nice pacific island and build my own winter retreat there.
37. Acquire an underground missile base from the US and make it nice.
38. Acquire some weapons and explosives for fun and distribute them to my retreats all around the world.
39. Stockpile the retreats with supplies and build bunkers.
40. Visit Antarctica and Cape Horn area, Greenland, Siberia, India, China, New Zealand, Jordan, Africa...
41. Build an internationally successful company that manufactures useful hardware, solves a number of practical problems and employs a number of good people.
42. Make a game that tells an interesting story, teaches a few things and presents practical problems requiring decent problem solving skills.
43. Work for a while as national park ranger, chemist, philosopher, politician and as a neuroscientist.
44. Work for a while on genetic research and related biological issues.
45. Work for a while for NASA, nvidia, AI research, weapons research, some startups...
46. Live a few years in the wilderness and a monestary.
47. Try living in different cities around the world.
48. Get a cabin from Norway from some nice remote mountain slope next to a small lake with a small airstrip where I can fly from Finland with my own small plane.
49. Get a nice place from somewhere along the Californian cost, not too far from San Francisco, but in a not too densely populated location either.
50. Write a book that will stir up people and teach them some critical thinking.

* The word "prove" used in colloquial sense here.

Though I didn't like this particular turn of events in that particular episode as in my opinion while it is not possible to avoid all mistakes given the initial conditions, I never the less believe that losing is due to mistakes that under perfect knowledge and logic should not exist. Choose your mistakes wisely and you'll be fine.

Unexamined life is not worth living is often taken to imply that rather than not being worth living some lives are simply missing something important. I would, however, argue that all lives are missing something important and on the other than that no single thing in any life is something so important that one should consider it essential. Most life has intrinsic value that is independent of its fullness.

It is often also said that the saying refers to one being able to say why things are important, right or wrong rather than simply insisting they are. However, this too is a bit arbitrary as in the end there never seems to be any why. We are what we are and want what we want, because that's the way we have come to be. Some goals are better served by some actions, but no action is right or wrong without reference to motivations which in end are arbitrary. Without such motivations we would not be humans, but rather simple machines simply waiting to do what we are told.

It is of course worthwhile to know how to influence other people by justifying actions and inactions for the betterment of humanity in general.

Origin of the power of quantum computing

Predicting the behavior of many particle quantum system requires solving the Schrödinger equation for a wavefunction that has complex dimensionality proportinal to the number of particles. Such a wavefunctions contains a value for all possible combinations of particle positions, i.e. it describes a kind of "multiverse" where everything conceivable happens simultaneously.

The system evolves as a whole and thus in order to predict phenomena such as entanglement, it is necessary to solve the time evolution of the "multiverse" rather than time evolution of separate wavefunctions corresponding to single particles. Single separate wavefunctions cannot yield predictions corresponding to observations in every circumstance.

The power of quantum computing is based on the complexity of the time evolution of this kind of multiparticle wavefunction. The classical computational resources required to do full simulations scale extremely poorly, but it is possible to demonstrate the idea with a simple code (below).

Full simulation of the time evolution of three electrons starting at rest relative to each other. Variable containing the wavefunction is a 6D grid of 26 complex points (3D double complex 26x26) and requires about 2.3 GBs of memory. Adding one more particle would make it 8D and require about 1556 GBs of memory. Simulating 30 particles in this configuration would require more bytes of memory than there are atoms in the observable universe. Many problems can be approximated effectively without conducting such simulations, but there are a few cases like quantum computers which cannot. Qubit are typically in a simple static grid rather than free space such as electrons in this simulation so larger number of qubits can still be simulated. The point at which quantum computers are expected to become more powerful than classical supercomputers (so called quantum supremacy) is a the moment 56 qubits.

n = 0.05;
N = 26;
x10 = -0.4; y10 =  0.4; vx1 = 0; vy1 = 0;
x20 =  0.4; y20 =  0.4; vx2 = 0; vy2 = 0;
x30 =  0.0; y30 = -0.4; vx3 = 0; vy3 = 0;
frequency_spread = 40;
psi = zeros(N, N, N, N, N, N, 'single');
V = zeros(N, N, N, N, N, N, 'single');
S = linspace(-1.0, 1.0, N);
for y1 = 1:N
    for x1 = 1:N
        for y2 = 1:N
            for x2 = 1:N
                for y3 = 1:N
                    for x3 = 1:N
                        envelope = exp(-frequency_spread*((S(x1)+x10)^2 + (S(y1)+y10)^2 + ...
                                                          (S(x2)+x20)^2 + (S(y2)+y20)^2 + ...
                                                          (S(x3)+x30)^2 + (S(y3)+y30)^2));
                        velocity = exp(-j*(vx1*S(x1) + vy1*S(y1) + ...
                                           vx2*S(x2) + vy2*S(y2) + ...
                                           vx3*S(x3) + vy3*S(y3)));
                        psi(x1, y1, x2, y2, x3, y3) =  envelope.*velocity;
                        V(x1, y1, x2, y2, x3, y3) = 1/sqrt(n + (S(x1)-S(x2))^2 + ...
                                                               (S(y1)-S(y2))^2)^2 + ...
                                                    1/sqrt(n + (S(x1)-S(x3))^2 + ...
                                                               (S(y1)-S(y3))^2)^2 + ...
                                                    1/sqrt(n + (S(x2)-S(x3))^2 + ...
                                                               (S(y2)-S(y3))^2)^2;
                    end
                end
            end
        end
    end
end
for t = 1:500
    PSI1 = squeeze(mean(mean(mean(mean(psi, 1), 2), 3), 4));
    PSI2 = squeeze(mean(mean(mean(mean(psi, 1), 2), 5), 6));
    PSI3 = squeeze(mean(mean(mean(mean(psi, 3), 4), 5), 6));
    a = abs(PSI1).^2; a = a/max(max(a)); b = abs(PSI2).^2; b = b/max(max(b));
    c = abs(PSI3).^2; c = c/max(max(c)); val = a + b + c;
    imwrite(val/max(max(val)), ['png/' num2str(t, '%04.f') '.png']);
    
    z = zeros(N, N, N, N, N, N, 'single');
    z(2:N-1, :, :, :, :, :) = z(2:N-1, :, :, :, :, :) + diff(psi, 2, 1);
    z(:, 2:N-1, :, :, :, :) = z(:, 2:N-1, :, :, :, :) + diff(psi, 2, 2);
    z(:, :, 2:N-1, :, :, :) = z(:, :, 2:N-1, :, :, :) + diff(psi, 2, 3);
    z(:, :, :, 2:N-1, :, :) = z(:, :, :, 2:N-1, :, :) + diff(psi, 2, 4);
    z(:, :, :, :, 2:N-1, :) = z(:, :, :, :, 2:N-1, :) + diff(psi, 2, 5);
    z(:, :, :, :, :, 2:N-1) = z(:, :, :, :, :, 2:N-1) + diff(psi, 2, 6);
    % that's just del2(psi)

    psi = psi + 0.005*j*(z - V.*psi);
end

--

mergesort(int a[], int n) { // O(n log n)
  int *b, *c;
  int x=0, y=0, z;
  if(n>3) {
    b = (int*)malloc(4*n/2);
    c = (int*)malloc(4*n/2);
    for(z=0; z<n/2; z++) {
      b[z] = a[z];
      c[z] = a[z+n/2];
    }
    mergesort(b, n/2);
    mergesort(c, n/2);
    for(z=0; z<n; z++)
      if(y<n/2 && x<n/2)
        if(b[x]>c[y]) a[z] = c[y++];
        else a[z] = b[x++];
      else
        if(y>x) a[z] = b[x++];
        else a[z] = c[y++];
    free(b);
    free(c);
  } else
    if(a[0]>a[1]) { // swap a[0] and a[1] (xor swap)                            
      a[0] = a[0]^a[1];
      a[1] = a[1]^a[0];
      a[0] = a[0]^a[1];
    }
}

--

Standard bicubic upscaler (4x).

A = imread('smb.png');
A = double(A)/255;

s = size(A);
S = 2;
for y = 1:S*s(2)-1
    for x = 1:S*s(1)-1
        r(x, y) = A(1+floor(x/S), 1+floor(y/S), 1);
        g(x, y) = A(1+floor(x/S), 1+floor(y/S), 2);
        b(x, y) = A(1+floor(x/S), 1+floor(y/S), 3);
    end
end

r = r + del2(r);
g = g + del2(g);
b = b + del2(b);

...repeat...

[x y] = gradient(r); r = r + abs(x) + abs(y);
[x y] = gradient(g); g = g + abs(x) + abs(y);
[x y] = gradient(b); b = b + abs(x) + abs(y);
[x y] = gradient(r); r = r - abs(x) - abs(y);
[x y] = gradient(g); g = g - abs(x) - abs(y);
[x y] = gradient(b); b = b - abs(x) - abs(y);

Nabla squared filter applied in two stages to nearest neighbourhood doubler and sharpen in post prosessing. Pixel edges in 45 degree angle are totally eliminated and regular raster becomes solid color.

Schrödinger sims, raytracing and Julia sets

Particle colliding with an obstacle.

N = 256;

[x y] = meshgrid(linspace(-1, 1, N), linspace(-1, 1, N));
psi = exp(-10*((x+0.5).^2 + (y+0.5).^2)).*exp(-j*(-50*x -50*y));
psi = psi/norm(psi);
V = zeros(N, N);
V(N/2-20:N/2-10, N/2:N/2+20) = 1;
V(N/2+10:N/2+20, N/2-10:N/2+10) = 1;
S = conv2(diff([0 0 1 0 0], 2), eye(N)); % del2 matrix form
for n = 1:N % H = V - del2
    Ex{n} = expm(-j*(diag(V(n, :)) - S(:, 2:end-1)));
    Ey{n} = expm(-j*(diag(V(:, n)) - S(:, 2:end-1)));
end
for t = 0:150 % time evolution psi(t+1) = psi(t)*exp(-j*H)
    for n = 1:N
        psi(:, n) = Ex{n}*psi(:, n);
    end
    for n = 1:N
        psi(n, :) = psi(n, :)*Ey{n};
    end
    A = abs(psi).^2/4.5e-4;
    cla; imagesc(A); caxis([0 1]); truesize; drawnow;
    imwrite(imresize(A, 0.5), ['png/' num2str(t+1, '%04.f') '.png']);
    % ffmpeg -r 30 -f image2 -i %04d.png -vcodec libx264 -crf 25 -pix_fmt yuv420p test.mp4
    % convert -resize 50% -layers Optimize -delay 3 *.png test.gif
end

Three particles scattering, simplified approximate solution only.

N = 64;
[x y] = meshgrid(linspace(-1, 1, N), linspace(-1, 1, N));
A = 1.4; B = -0.2; C = 1.2;
psi{1} = exp(-50*((x-0.5*A-B).^2 + (y-0.4*A).^2)).*exp(-j*( 20*x*C + 10*y*C));
psi{2} = exp(-50*((x+0.2*A-B).^2 + (y+0.3*A).^2)).*exp(-j*(-20*x*C - 10*y*C));
psi{3} = exp(-50*((x-0.4*A-B).^2 + (y+0.3*A).^2)).*exp(-j*( 10*x*C - 20*y*C));
v = fft2(0.3./(x.^2+y.^2+0.01));
for x = 1:3
    psi{x} = psi{x}/norm(psi{x});
end
in = 1;
for t = 0:15000
    if mod(t, 150)==0
        B = 0;
        for x = 1:3
            A = abs(psi{x}).^2;
            B = B + A/max(max(A));
        end
        figure(1); cla; imagesc(imresize(B, 4)); truesize; drawnow;
        imwrite(B, ['png/' num2str(in, '%04.f') '.png']);
        in = in + 1
    end
    for x = 1:3 % fast 2D convolution with fft
        V{x} = fftshift(ifft2(v.*fft2(conj(psi{x}).*psi{x})));
    end
    for x = 1:3 % separate wave functions
        W = zeros(N, N);
        for y = 1:3
            if y~=x
                W = W + V{y};
            end
        end
        psi{x} = psi{x} + 1e-3*j*(4*del2(psi{x}) - W.*psi{x});
    end   
end

Traces ever pixel of the screen to z-direction until ray within the object and
then chooses color based on surface normal and distance.
The edges to infinity are funny due to "cutting corners" in the way the normal is computed.

/* gcc test.c -lGL -lSDL */
#include <GL/gl.h>
#include <SDL/SDL.h>

char *d =
  "float t = gl_Color.x*1000.;"
  "vec3 scene(vec2 p) {"
  "  float z = 5.;"
  "  mat3 X = mat3(1., 0., 0., 0., cos(2.*t), sin(2.*t), 0., -sin(2.*t), cos(2.*t));"
  "  mat3 Y = mat3(cos(3.*t), 0., -sin(3.*t), 0., 1., 0., sin(3.*t), 0., cos(3.*t));"
  "  vec3 c = X*Y*vec3(1.5*cos(.5*5.*t), 1.5*sin(.5*7.*t), 7.);"
  "  while(z<10.) {"
  "    vec3 r = X*Y*(vec3(sqrt(z)*p.x, sqrt(z)*p.y, z));"
  "    if(abs(r.x-c.x)<.8 && abs(r.y-c.y)<.8 && abs(r.z-c.z)<.8) return r;"
  "    if(distance(r, c)<1.) return r;"
  "    z += .05;"
  "  }"
  "  return vec3(10.);"
  "}"
  "void main() {"
  "  gl_FragColor = vec4(0.);"
  "  vec2 p = vec2(gl_FragCoord.x/256.-1., gl_FragCoord.y/256.-1.);"
  "  vec3 q = scene(p);"
  "  if(length(q)<10.) {"
  "    vec3 dx = scene(vec2(p.x + .05, p.y)) - scene(vec2(p.x - .05, p.y));"
  "    vec3 dy = scene(vec2(p.x, p.y + .05)) - scene(vec2(p.x, p.y - .05));"
  "    float a = 1.5*dot(normalize(cross(dx, dy)), normalize(q))/sqrt(length(q));"
  "    gl_FragColor = vec4(a);"
  "  }"
  "}";

main(c) {
  SDL_Event e;
  SDL_SetVideoMode(512, 512, 0, SDL_OPENGL);
  glUseProgram(glCreateShaderProgramv(GL_FRAGMENT_SHADER, 1, &d));

/*
  c = 0;
  GLubyte **pixels = malloc(3*512*512);
  int cur;
  char filename[80];

*/

  while(e.type!=SDL_QUIT) {
    glColor3us(c++, 0, 0);
    glRecti(-1, -1, 1, 1);
    SDL_GL_SwapBuffers();
    SDL_PollEvent(&e);
/*
    sprintf(filename, "%04d.ppm", c);
    FILE *f = fopen(filename, "w");
    fprintf(f, "P3\n%d %d\n%d\n", 512, 512, 255);
    *pixels = realloc(*pixels, 3*512*512);
    glReadPixels(0, 0, 512, 512, GL_RGB, GL_UNSIGNED_BYTE, *pixels);
    for(int i=0; i<512; i++) {
        for(int j=0; j<512; j++) {
            cur = 3*((512-i-1)*512+j);
            fprintf(f, "%3d %3d %3d ", (*pixels)[cur], (*pixels)[cur+1], (*pixels)[cur+2]);
        }
        fprintf(f, "\n");

    }

*/

  }
}

Zoom into the Julia set (one of them).

N = 256; B = 0.75;

X = linspace(-1, 1, N);
Y = linspace(-1, 1, N);
while true
    B = B*sqrt(sqrt(0.75))
    y0 = 0.636;
    x0 = 0.107;
    X = linspace(x0-B, x0+B, N);
    Y = linspace(y0-B, y0+B, N);

    for y = 1:N

        for x = 1:N

            % z = 0.0;                  % Mandelbrot

            % z0 = X(x) + Y(y)*j;       % Mandelbrot

            z = X(x) + Y(y)*j;          % Julia

            z0 = 0.5 + 0.5*j;           % Julia

            iteration = 0;

            while abs(z) < 2 && iteration < 128

                z = z^2 + z0;

                iteration = iteration + 1;

            end

            s(x, y) = log(iteration);

        end

    end

    imagesc(s); drawnow;
end

#include <stdio.h>
#include <math.h>
main() {
  char filename[80];
  int N = 128, x, y, c=1;
  float f[128*128];  float p[128*128];
  float u[128*128];  float v[128*128];
  float ux[128*128]; float vx[128*128];
  float uy[128*128]; float vy[128*128];
  float u2[128*128]; float v2[128*128];
  float us[128*128]; float vs[128*128];
  float px[128*128]; float py[128*128];
  for(x=0; x<128*128; x++) {
    f[x] = 0.0;    p[x] = 0.0;
    u[x] = 0.0;    v[x] = 0.0;
    ux[x] = 0.0;   vx[x] = 0.0;
    uy[x] = 0.0;   vy[x] = 0.0;
    u2[x] = 0.0;   v2[x] = 0.0;
    us[x] = 0.0;   vs[x] = 0.0;
    px[x] = 0.0;   py[x] = 0.0;  }
  for(int t=0; t<4500; t++) {
    for(y=-5; y<5; y++)
      for(x=-5; x<5; x++) {
        u[(int)(0.75*N+x)+(N/2+y)*N] = 0.15;
        v[(int)(0.75*N+x)+(N/2+y)*N] = 0.0;
      }
    for(y=1; y<N-1; y++)
      for(x=1; x<N-1; x++) {
        ux[x+y*N] = 0.50*(u[(x+1)+y*N] -                u[(x-1)+y*N]);
        vx[x+y*N] = 0.50*(v[(x+1)+y*N] -                v[(x-1)+y*N]);
        uy[x+y*N] = 0.50*(u[x+(y+1)*N] -                u[x+(y-1)*N]);
        vy[x+y*N] = 0.50*(v[x+(y+1)*N] -                v[x+(y-1)*N]);
        u2[x+y*N] = 0.25*(u[(x+1)+y*N] - 2.0*u[x+y*N] + u[(x-1)+y*N]);
        v2[x+y*N] = 0.25*(v[(x+1)+y*N] - 2.0*v[x+y*N] + v[(x-1)+y*N]);
        u2[x+y*N]+= 0.25*(u[x+(y+1)*N] - 2.0*u[x+y*N] + u[x+(y-1)*N]);
        v2[x+y*N]+= 0.25*(v[x+(y+1)*N] - 2.0*v[x+y*N] + v[x+(y-1)*N]);
      }
    for(x=0; x<N*N; x++) {
      us[x] = u[x] + u[x]*ux[x] + v[x]*uy[x] + u2[x]/20.0 - px[x]/20.0 + 0.00000;
      vs[x] = v[x] + u[x]*vx[x] + v[x]*vy[x] + v2[x]/20.0 - py[x]/20.0 - 0.00005;
    }
    for(y=1; y<N-1; y++)
      for(x=1; x<N-1; x++) {
        f[x+y*N] = 0.5*(us[(x+1)+y*N] - us[(x-1)+y*N] + vs[x+(y+1)*N] - vs[x+(y-1)*N]);
        p[x+y*N] = 0.25*(p[(x+1)+y*N] + p[x+(y+1)*N] + p[(x-1)+y*N] + p[x+(y-1)*N] - f[x+y*N]);
    }
    for(y=1; y<N-1; y++)
      for(x=1; x<N-1; x++) {
        px[x+y*N] = 0.5*(p[(x+1)+y*N] - p[(x-1)+y*N]);
        py[x+y*N] = 0.5*(p[x+(y+1)*N] - p[x+(y-1)*N]);
         u[x+y*N] = us[x+y*N] - px[x+y*N];
         v[x+y*N] = vs[x+y*N] - py[x+y*N];
      }
    if(t%40==0) {
      printf("%d\n", t);
      sprintf(filename, "ppm/%04d.pgm", c);
      FILE *f = fopen(filename, "w");
      fprintf(f, "P2\n%d %d\n%d\n", 128, 128, 255);
      c++;
      for(int j=0; j<128; j++) {
          for(int i=0; i<128; i++) {
             float val = 1500.0*sqrt( u[i+j*N]*u[i+j*N] + v[i+j*N]*v[i+j*N]);
             if(val>255.0) val=255.0;
             fprintf(f, "%3d ", (unsigned char)val);
          }
          fprintf(f, "\n");
      }
      close(f);
    }
  }
}

KISS

-Safety, speed, usability, etc. all generally speaking originate from the same source, the system doing as little as possible in the simplest and understandable way possible while achieving the goal of doing a lot of the same thing it's good at.

-Transparency: the system must be transparent and open. It must not do anything behind the users back. The user must be able to view all that is going on in their system and understand it.

-Non-redundancy: in UI there ought to be only one way of doing things so as not to confuse the user or waste their time by making them learn many ways of doing the same thing.

-Limited abstraction: working close to the hardware while maintaining a level of readability and functionality. Unnecessary abstraction leads to waste of time, less speed, worse security and confusion.

-Stay out of the way of the user, never popup anything unless requested, never run anything in the background unless specifically requested.

Visit Russia, find subs, observe Militsiya, get a cold (after two years of health).

WM design ideas:

Applications are allowed to bring themselves on top only once (when they are started within certain time). After that the user must ALT+TAB to bring them on top.

OS/Compiler stuff:

1) Compile a static busybox and copy the files to a HDD image (after making the partition):
mount -o offset=32256 sda.bin tmp
2) mknod a few devices
3) Compile linux kernel

Tools these days...

4) qemu-system-x86_64 -m 512 -kernel bzImage -append "root=/dev/sda1" -enable-kvm -hda sda.bin

5) Modify tcc so it can be compiled into static binary (./configure --extra-cflags=-static), improve the assembler by adding syscall, use statically linked uClibc for startup code (plan: write the startup code entirely from scratch).

x86_64-asm.h: DEF_ASM_OP0(syscall, 0x0f05)
tccpp.c: get rid of ldexp and make your own d = d*((1<<((int)(exp_val-frac_bits))));
tccelf.c: get rid of all things related to libtcc1.a and crt?.o
tcc.c: add missing (dummy) function void __assert_fail (const char *assertion, const char *file, unsigned int line, const char *function) { abort(); }

tcc -DCONFIG_TCC_STATIC -o tcc.o -c tcc.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o libtcc.o -c libtcc.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o tccpp.o -c tccpp.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o tccgen.o -c tccgen.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o tccelf.o -c tccelf.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o tccasm.o -c tccasm.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o tccrun.o -c tccrun.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o x86_64-gen.o -c x86_64-gen.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o i386-asm.o -c i386-asm.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o libtcc1.o -c libtcc1.c -DTCC_TARGET_X86_64 -I./include -static
tcc -DCONFIG_TCC_STATIC -o alloca86_64.o -c alloca86_64.S -DTCC_TARGET_X86_64 -I./include -static

gcc -nostdlib -static -o tcc libtcc1.o alloca86_64.o tcc.o libtcc.o tccpp.o tccgen.o tccelf.o tccasm.o tccrun.o x86_64-gen.o i386-asm.o ../uClibc-0.9.33/lib/crt1.o ../uClibc-0.9.33/lib/crti.o ../uClibc-0.9.33/lib/crtn.o ../uClibc-0.9.33/lib/libc.a
strip tcc

6) End up with standalone 290 kB static C-compiler/assembler.

Boots in a second.

void* syscall(void* syscall_number, void* param1, void* param2, void* param3) {

  __asm__ __volatile__(
    "mov %0, %%rax\n"
    "mov %1, %%rdi\n"
    "mov %2, %%rsi\n"
    "mov %3, %%rdx\n"
    "syscall\n"
    :
    : "g" (syscall_number), "g" (param1), "g" (param2), "g" (param3)
    : "rax", "rdi", "rsi", "rdx");
}

typedef unsigned long int uintptr;
typedef long int intptr;

static intptr write(int fd, void const* data, uintptr nbytes) {
  syscall((void*)1, (void*)(intptr)fd, (void*)data, (void*)nbytes);
}

int main(int argc, char* argv[]) {
    write(1, "Hello, World!\n", 15);

    return 0;
}

_start() {
  __asm__ __volatile__(
    "call main\n"
    "mov $60, %rax\n"
    "mov 0, %rdi\n"
    "syscall\n");
}

7) Compile hello world on a few MB linux installation on a 64 bit CPU.

Belonging, purpose, transcendence and storytelling...

Navier-Stokes, good luck making the code shorter...

N  = 192; dt = 0.25; p = zeros(N, N);
u = zeros(N, N); v = zeros(N, N); ustar = zeros(N, N); vstar = zeros(N, N);
k = pi*[0:(N/2-1) (-N/2):(-1)]; [x y]  = meshgrid(k, k); ds = -(x.^2 + y.^2);
ds(1,1) = 1; px = 0.0; py = 0.0; gx = 0.0; gy = 0.0;
for t = 0:1000
    % constant velocity zone
    u(N/2-10:N/2+10, N/2-10:N/2+10) = 0.2;

    % solve intermediate velocity U*
    [ux uy] = gradient(u);
    [vx vy] = gradient(v);
    ustar = u + (u.*ux + v.*uy + del2(u)/20 - px/20 + gy)*dt;
    vstar = v + (u.*vx + v.*vy + del2(v)/20 - py/20 + gx)*dt;

    % solve Poisson equation: del2(p) = f
    fhat = fft2(divergence(ustar, vstar));
    p = N^2*real(ifft2(fhat./ds));

    % enforce del(U) = 0
    [px py] = gradient(p);
    u = ustar - px*dt;
    v = vstar - py*dt;

    cla; surface(sqrt(u.^2 + v.^2), 'edgecolor', 'none'); drawnow;
end

Espoo is full of furry animals.

Just a basic wave equation...

u  = zeros(192, 192);
du = zeros(192, 192);
[x y] = meshgrid(linspace(-1, 1, 192), linspace(-1, 1, 192));
for t = 0:0.002:0.4
    u = u + exp(-1000*((x-0.7*sin(7*t)).^2+(y-0.7*sin(11*t)-0.2).^2));
    u = u + 0.5*exp(-1000*((x+0.7*sin(5*t)).^2+(y+0.7).^2))*sin(3*31*t);

    du = du + del2(u);   
    u = u + du;
       
    cla; imagesc(u); drawnow;
end

...