Skip to content

Instantly share code, notes, and snippets.

@skeeto skeeto/Makefile

Created Sep 19, 2017
Embed
What would you like to do?
Carpet Genetic Algorithm
#include <float.h>
#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <string.h>
#define TARGET 81
#define COLORS 32
#define RULE(r, n, p) (0x1f & (r[n] >> ((p) * 5)))
#define DX(i) ((int)((0x0489a621UL >> (4 * (i) + 0)) & 3) - 1)
#define DY(i) ((int)((0x0489a621UL >> (4 * (i) + 2)) & 3) - 1)
#define COUNTOF(a) (sizeof(a) / sizeof(*a))
static uint32_t
pcg32(uint64_t *s)
{
*s = *s * UINT64_C(0x9b60933458e17d7d) + UINT64_C(0xd737232eeccdf7ed);
int shift = 29 - (*s >> 61);
return *s >> shift;
}
static int
pgm_load(unsigned char *buf)
{
int w, h, d;
if (scanf("P2 %d %d %d", &w, &h, &d) != 3)
return 0;
if (w != TARGET || h != TARGET || d != 31)
return 0;
for (int i = 0; i < TARGET * TARGET; i++) {
int v;
scanf("%d", &v);
if (v < 0 || v > d)
return 0;
buf[i] = v;
}
return 1;
}
static void
pgm_write(unsigned char *buf, int scale, FILE *o)
{
fprintf(o, "P2\n%d %d\n%d\n", TARGET * scale, TARGET * scale, COLORS);
for (int y = 0; y < TARGET * scale; y++)
for (int x = 0; x < TARGET * scale; x++)
fprintf(o, "%d\n", buf[(y / scale) * TARGET + (x / scale)]);
}
static void
ppm_write(unsigned char *buf, int scale, FILE *o)
{
fprintf(o, "P6\n%d %d\n255\n", TARGET * scale, TARGET * scale);
for (int y = 0; y < TARGET * scale; y++) {
for (int x = 0; x < TARGET * scale; x++) {
int v = buf[(y / scale) * TARGET + (x / scale)];
v = v * 255 / (COLORS - 1);
fputc(v, o);
fputc(v, o);
fputc(v, o);
}
}
}
static double
score(unsigned char *a, unsigned char *b)
{
double error = 0.0;
for (int y = 0; y < TARGET; y++) {
for (int x = 0; x < TARGET; x++) {
int c = 1;
int amean = a[y * TARGET + x];
int bmean = b[y * TARGET + x];
for (int i = 0; i < 8; i++) {
int xx = x + DX(i);
int yy = y + DY(i);
if (xx >= 0 && xx < TARGET && yy >= 0 && yy < TARGET) {
amean += a[yy * TARGET + xx];
bmean += b[yy * TARGET + xx];
c++;
}
}
double delta = amean / (double)c - bmean / (double)c;
error += delta * delta;
}
}
return error;
}
static void
render(unsigned char *buf, unsigned long long *rules)
{
for (int y = 0; y < TARGET; y++) {
for (int x = 0; x < TARGET; x++) {
int c = 0;
int px = x;
int py = y;
int div = TARGET / 3;
for (int i = 1; div > 0; i++) {
c = RULE(rules, c, py / div * 3 + px / div);
px %= div;
py %= div;
div /= 3;
}
buf[y * TARGET + x] = c;
}
}
}
static void
random_rules(unsigned long long *rules, uint64_t *rng)
{
for (int i = 0; i < COLORS; i++) {
uint64_t high = pcg32(rng);
uint64_t low = pcg32(rng);
rules[i] = (high << 32) | low;
}
}
static void
breed(unsigned long long *o,
unsigned long long *a,
unsigned long long *b,
uint64_t *rng)
{
uint32_t select = pcg32(rng);
for (int i = 0; i < COLORS; i++)
o[i] = (select >> i) & 1 ? a[i] : b[i];
int mutations = pcg32(rng) % 32;
for (int i = 0; i < mutations; i++) {
select = pcg32(rng);
int r = (select >> 0) % COLORS;
int c = (select >> 5) % 9;
int t = (select >> 9) % COLORS;
uint64_t mask = UINT64_C(0x1f) << (c * 5);
o[r] = (o[r] & ~mask) | ((uint64_t)t << (c * 5));
}
}
static int
dblcmp(const void *pa, const void *pb)
{
double a = *(double *)pa;
double b = *(double *)pb;
if (a < b)
return -1;
else if (a > b)
return 1;
return 0;
}
static void
rule_write(unsigned long long *rule, FILE *o)
{
fputs("static const unsigned long long rules[] = {\n", o);
for (int i = 0; i < COLORS; i++) {
if (i % 4 == 0)
fputs(" ", o);
unsigned long long v = rule[i] & 0x1fffffffffff;
fprintf(o, "0x%012llx,%c", v, " \n"[i % 4 == 3]);
}
fputs("};\n", o);
for (int i = 0; i < COLORS; i++) {
fprintf(o, "// %d ->", i);
for (int c = 0; c < 9; c++) {
int p = RULE(rule, i, c);
fprintf(o, " %d%c", p, ",\n"[c == 8]);
}
}
}
int main(void)
{
uint64_t rng[] = {0x9e8480dd162324e1};
unsigned long long rules[200][COLORS];
unsigned char input[TARGET * TARGET];
if (!pgm_load(input)) {
fputs("Invalid input\n", stderr);
return 1;
}
for (size_t i = 0; i < COUNTOF(rules); i++)
random_rules(rules[i], rng);
FILE *video = fopen("video.ppm", "wb");
double global_best = DBL_MAX;
for (long long g = 0; ; g++) {
struct {
double score;
size_t i;
} scores[COUNTOF(rules)];
unsigned long long best[10][COLORS];
/* Find the 10 best rulesets */
for (size_t i = 0; i < COUNTOF(rules); i++) {
unsigned char buf[TARGET * TARGET];
render(buf, rules[i]);
scores[i].score = score(input, buf);
scores[i].i = i;
}
qsort(scores, COUNTOF(scores), sizeof(*scores), dblcmp);
/* Breed next generation */
for (int i = 0; i < COUNTOF(best); i++)
memcpy(best[i], rules[scores[i].i], sizeof(best[i]));
for (int i = 0; i < COUNTOF(rules) - 5; i++) {
uint32_t select = pcg32(rng);
int a = (select >> 0) % 10;
int b = (select >> 16) % 10;
breed(rules[i], best[a], best[b], rng);
}
for (int i = COUNTOF(rules) - 5; i < COUNTOF(rules); i++)
random_rules(rules[i], rng);
/* Report on progress */
if (scores[0].score < global_best) {
global_best = scores[0].score;
/* Write out the current best image */
unsigned char buf[TARGET * TARGET];
FILE *image = fopen("_progress.pgm", "w");
render(buf, best[0]);
pgm_write(buf, 10, image);
fclose(image);
rename("_progress.pgm", "progress.pgm");
/* Write out best image as video frame */
ppm_write(buf, 10, video);
fflush(video);
/* Write out the ruleset */
FILE *save = fopen("best.c", "w");
rule_write(best[0], save);
fclose(save);
}
printf("%12lld %f %f\n", g, global_best, scores[0].score);
}
}
CC = c99
CFLAGS = -Wall -Wextra -Ofast
all: genetic render
genetic: genetic.c
$(CC) $(CFLAGS) $(LDFLAGS) -o $@ genetic.c $(LDLIBS)
render: render.c
$(CC) $(CFLAGS) $(LDFLAGS) -o $@ render.c $(LDLIBS)
#include <stdio.h>
#define TARGET 81
#define SCALE 10
#define RULE(n, p) (0x1f & (rules[n] >> ((p) * 5)))
static const unsigned long long rules[] = {
0x12cff8e3eefd, 0x0240c6930008, 0x1ff73beef75b, 0x16dfefffebfd,
0x1e6ee3c67375, 0x1d2b8e100422, 0x1d17fbee6b0a, 0x1aed8664928b,
0x0e4918f5fc67, 0x0509ec4e0c41, 0x0167184a847f, 0x02824300785a,
0x0002a000a2e7, 0x17ef6e67ce76, 0x1af987ff3def, 0x1af8fbef94a5,
0x14281121307f, 0x0cb535de8846, 0x1feffa65987f, 0x002c400099a8,
0x0508bab67842, 0x1dd9b7457416, 0x080801111522, 0x02875e3f8bff,
0x02f17a988ff0, 0x0e7ba2f8c451, 0x1dff79ff739f, 0x1f7ff0e6fe3d,
0x07e995a15160, 0x1dfffbfffffd, 0x1be73fee7bff, 0x1dff7daff7da,
};
int main(void)
{
printf("P2\n%d %d\n29\n", TARGET * SCALE, TARGET * SCALE);
for (int y = 0; y < TARGET * SCALE; y++) {
for (int x = 0; x < TARGET * SCALE; x++) {
int c = 0;
int px = x / SCALE;
int py = y / SCALE;
int div = TARGET / 3;
for (int i = 1; div > 0; i++) {
c = RULE(c, py / div * 3 + px / div);
px %= div;
py %= div;
div /= 3;
}
printf("%d\n", c);
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.