Files
LithosAnanake/src/compudynamics.c
T

383 lines
14 KiB
C
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
/*
StarForth — Steady-State Virtual Machine Runtime
Copyright (c) 20232025 Robert A. James
All rights reserved.
This file is part of the StarForth project.
Licensed under the StarForth License, Version 1.0 (the "License");
you may not use this file except in compliance with the License.
You may obtain a copy of the License at:
https://github.com/star.4th@proton.me/StarForth/LICENSE.txt
This software is provided "AS IS", WITHOUT WARRANTY OF ANY KIND,
express or implied, including but not limited to the warranties of
merchantability, fitness for a particular purpose, and noninfringement.
See the License for the specific language governing permissions and
limitations under the License.
*/
#include "compudynamics.h"
#include <string.h>
#include <stdlib.h>
/* ============================================================================
* Per-Level Tuning Accessors
* ============================================================================
*
* Returned by value, not exposed as addressable globals -- see the
* "Tuning knobs are returned by value" comment in compudynamics.h for why.
*/
/* Matches the constants ssm_jacquard.c/.h used before this module existed
* (SSM_SCORE_MAX, SSM_SCORE_DECAY_FACTOR_Q16, SSM_UCB_K, SSM_REWARD_GAIN,
* SSM_JOINT_WEIGHT_Q16, SSM_ANOVA_WEIGHT_Q16, SSM_ENTROPY/CV/TEMPORAL_DECAY
* thresholds) -- a lift, not a retune. id_window_depth=32 and
* locality_lookback=8 are new: they size the entity-ID history the
* deterministic locality signal replaces temporal_decay with. */
CDTuning cd_tuning_word(void)
{
CDTuning t;
/* classifier */
t.id_window_depth = 32u;
t.locality_lookback = 8u;
t.diversity_high_threshold_q16 = 49152u; /* 0.75 * 65536 */
t.volatility_high_threshold_q16 = 9830u; /* 0.15 * 65536 */
t.locality_high_threshold_q16 = 32768u; /* 0.5 * 65536 */
/* bandit */
t.num_configs = 128u; /* 7-bit L1-L7 */
t.num_regimes = 8u; /* 3-bit diversity/volatility/locality */
t.score_max = 65535u;
t.score_min = 0u;
t.score_decay_factor_q16 = 64881u; /* 0.99 * 65536 */
t.ucb_k = 6554u; /* 10% of score_max */
t.reward_gain = 16384u; /* 25% of score_max */
t.joint_weight_q16 = 49152u; /* 0.75 */
t.anova_weight_q16 = 16384u; /* 0.25 */
return t;
}
/* VM level: classifier-side parameters only, matching capsule_vm_physics.c's
* existing VM_FLEET_WINDOW_DEPTH so a future migration has a ready-made
* starting point. Bandit fields are zero and unused -- see compudynamics.h. */
CDTuning cd_tuning_vm(void)
{
CDTuning t;
/* classifier */
t.id_window_depth = 64u; /* matches VM_FLEET_WINDOW_DEPTH */
t.locality_lookback = 16u;
t.diversity_high_threshold_q16 = 49152u;
t.volatility_high_threshold_q16 = 9830u;
t.locality_high_threshold_q16 = 32768u;
/* bandit -- unused at VM level today */
t.num_configs = 0u;
t.num_regimes = 0u;
t.score_max = 0u;
t.score_min = 0u;
t.score_decay_factor_q16 = 0u;
t.ucb_k = 0u;
t.reward_gain = 0u;
t.joint_weight_q16 = 0u;
t.anova_weight_q16 = 0u;
return t;
}
/* ============================================================================
* Integer Helpers (isqrt, ln approximation)
* ============================================================================
*/
static uint32_t cd_isqrt32(uint32_t n)
{
uint32_t x, y;
if (n == 0u) return 0u;
x = n;
y = (x + 1u) / 2u;
while (y < x) {
x = y;
y = (x + n / x) / 2u;
}
return x;
}
static uint64_t cd_isqrt64(uint64_t n)
{
uint64_t x, y;
if (n == 0u) return 0u;
x = n;
y = (x + 1u) / 2u;
while (y < x) {
x = y;
y = (x + n / x) / 2u;
}
return x;
}
/* floor(log2(n)) * ln(2) * 65536, integer shifts only. See ssm_jacquard.c's
* former ln_approx_q16 -- identical algorithm, moved here so the generic
* bandit doesn't depend on ssm_jacquard.c. */
static uint32_t cd_ln_approx_q16(uint32_t n)
{
uint32_t bits, tmp;
if (n <= 1u) return 0u;
bits = 0u;
tmp = n >> 1u;
while (tmp > 0u) { tmp >>= 1u; bits++; }
return bits * 45426u;
}
/* ============================================================================
* Deterministic Regime Classification
* ============================================================================
*/
CDRegimeMetrics cd_classify_ids(const uint32_t *ids, uint32_t count,
uint32_t locality_lookback)
{
CDRegimeMetrics m;
uint32_t n, i, j, distinct, lookback;
uint8_t counted[CD_MAX_CLASSIFY_DEPTH];
uint32_t per_id_count[CD_MAX_CLASSIFY_DEPTH];
m.diversity_q16 = 0u;
m.volatility_q16 = 0u;
m.locality_q16 = 0u;
if (!ids || count == 0u) return m;
n = (count > CD_MAX_CLASSIFY_DEPTH) ? CD_MAX_CLASSIFY_DEPTH : count;
memset(counted, 0, n * sizeof(counted[0]));
/* --- Diversity: distinct IDs / n, and per-ID touch counts for the
* volatility pass below (single O(n^2) sweep, n <= 64) --- */
distinct = 0u;
for (i = 0u; i < n; i++) {
uint32_t c;
if (counted[i]) continue;
c = 1u;
counted[i] = 1u;
for (j = i + 1u; j < n; j++) {
if (ids[j] == ids[i]) { c++; counted[j] = 1u; }
}
per_id_count[distinct] = c;
distinct++;
}
m.diversity_q16 = (uint32_t)(((uint64_t)distinct * 65536u) / n);
/* --- Volatility: coefficient of variation of per-ID touch counts,
* computed in Q8 to avoid overflow, result in Q16 --- */
if (distinct > 0u) {
uint32_t mean_q8 = (n << 8u) / distinct;
uint64_t sum_var_q16 = 0u;
uint64_t stddev_q8, cv_q16;
for (i = 0u; i < distinct; i++) {
int64_t diff_q8 = (int64_t)(per_id_count[i] << 8u) - (int64_t)mean_q8;
uint64_t diff2_q16 = (uint64_t)(diff_q8 * diff_q8);
sum_var_q16 += diff2_q16;
}
stddev_q8 = cd_isqrt64(sum_var_q16 / distinct);
cv_q16 = (mean_q8 > 0u) ? ((stddev_q8 * 65536u) / mean_q8) : 0u;
m.volatility_q16 = (cv_q16 > 0xFFFFFFFFu) ? 0xFFFFFFFFu : (uint32_t)cv_q16;
}
/* --- Locality: fraction of touches (from lookback..n-1) that repeat
* an ID seen within the preceding lookback touches --- */
lookback = (locality_lookback > n) ? n : locality_lookback;
if (n > lookback && lookback > 0u) {
uint32_t matches = 0u;
uint32_t denom = n - lookback;
for (i = lookback; i < n; i++) {
for (j = i - lookback; j < i; j++) {
if (ids[j] == ids[i]) { matches++; break; }
}
}
m.locality_q16 = (uint32_t)(((uint64_t)matches * 65536u) / denom);
}
return m;
}
/* ============================================================================
* Generic UCB1 Config-Space Bandit
* ============================================================================
*/
int cd_config_table_init(CDConfigTable *table, const CDTuning *tuning)
{
if (!table || !tuning) return -1;
memset(table, 0, sizeof(*table));
if (tuning->num_configs == 0u || tuning->num_regimes == 0u) return -1;
table->entries = (CDConfigEntry *)malloc(tuning->num_configs * sizeof(CDConfigEntry));
table->regime_scores = (uint32_t *)malloc(
(size_t)tuning->num_regimes * tuning->num_configs * sizeof(uint32_t));
table->regime_ticks = (uint32_t *)malloc(
(size_t)tuning->num_regimes * tuning->num_configs * sizeof(uint32_t));
table->total_regime_ticks = (uint32_t *)malloc(tuning->num_regimes * sizeof(uint32_t));
if (!table->entries || !table->regime_scores || !table->regime_ticks ||
!table->total_regime_ticks) {
cd_config_table_free(table);
return -1;
}
table->num_configs = tuning->num_configs;
table->num_regimes = tuning->num_regimes;
return 0;
}
void cd_config_table_free(CDConfigTable *table)
{
if (!table) return;
free(table->entries);
free(table->regime_scores);
free(table->regime_ticks);
free(table->total_regime_ticks);
memset(table, 0, sizeof(*table));
}
void cd_config_table_seed(CDConfigTable *table, const CDTuning *tuning,
const uint32_t *initial_scores, uint8_t initial_config)
{
uint32_t c, r;
if (!table || !table->entries || !tuning || !initial_scores) return;
for (c = 0u; c < table->num_configs; c++) {
CDConfigEntry *e = &table->entries[c];
e->config_bits = (uint8_t)c;
e->has_prev = 0u;
e->_pad[0] = 0u;
e->_pad[1] = 0u;
e->benefit_score = initial_scores[c];
e->tick_count = 0u;
e->prev_window = 0u;
e->prev_locality_q16 = 0u;
for (r = 0u; r < table->num_regimes; r++) {
table->regime_scores[r * table->num_configs + c] = initial_scores[c];
table->regime_ticks[r * table->num_configs + c] = 0u;
}
}
table->current_config = initial_config;
table->current_regime = 0u;
for (r = 0u; r < table->num_regimes; r++) {
table->total_regime_ticks[r] = 0u;
}
(void)tuning;
}
uint8_t cd_config_table_tick(CDConfigTable *table, const CDTuning *tuning,
uint8_t regime, uint32_t anova_stability_q16,
int use_joint_convergence,
uint32_t current_window, uint32_t current_locality_q16)
{
uint32_t r, c;
CDConfigEntry *entry;
uint32_t combined_stability;
uint32_t *score_cell;
uint32_t old_score, decayed, new_score;
int64_t delta_score, new_s;
uint32_t total, ln_total, best_c, best_ucb, ci;
if (!table || !table->entries || !tuning) return 0u;
table->current_regime = regime;
r = (uint32_t)regime;
c = (uint32_t)table->current_config;
entry = &table->entries[c];
if (use_joint_convergence && entry->has_prev) {
uint32_t dw_raw = (current_window > entry->prev_window)
? (current_window - entry->prev_window)
: (entry->prev_window - current_window);
uint32_t dl_raw = (current_locality_q16 > entry->prev_locality_q16)
? (current_locality_q16 - entry->prev_locality_q16)
: (entry->prev_locality_q16 - current_locality_q16);
uint64_t dw_q16 = (entry->prev_window > 0u)
? (((uint64_t)dw_raw * 65536u) / entry->prev_window)
: 0u;
uint64_t dl_q16 = dl_raw; /* already Q16, bounded to [0,65536] */
uint64_t sum_sq, joint_err;
uint32_t stability;
if (dw_q16 > 65536u) dw_q16 = 65536u;
if (dl_q16 > 65536u) dl_q16 = 65536u;
/* joint_error = sqrt(dw^2 + dl^2) / sqrt(2), normalised to [0,65536] */
sum_sq = dw_q16 * dw_q16 + dl_q16 * dl_q16;
joint_err = (cd_isqrt64(sum_sq) * 65536u) / 92682u; /* sqrt(2)*65536 ~= 92682 */
if (joint_err > 65536u) joint_err = 65536u;
stability = (uint32_t)(65536u - (uint32_t)joint_err);
combined_stability = (uint32_t)(
((uint64_t)tuning->joint_weight_q16 * stability +
(uint64_t)tuning->anova_weight_q16 * anova_stability_q16) / 65536u
);
} else {
combined_stability = anova_stability_q16;
}
delta_score = ((int64_t)combined_stability - 32768) *
(int64_t)tuning->reward_gain / 32768;
score_cell = &table->regime_scores[r * table->num_configs + c];
old_score = *score_cell;
decayed = (uint32_t)(((uint64_t)old_score * tuning->score_decay_factor_q16) / 65536u);
new_s = (int64_t)decayed + delta_score;
if (new_s < (int64_t)tuning->score_min) new_score = tuning->score_min;
else if (new_s > (int64_t)tuning->score_max) new_score = tuning->score_max;
else new_score = (uint32_t)new_s;
*score_cell = new_score;
table->regime_ticks[r * table->num_configs + c]++;
if (table->total_regime_ticks[r] < 0xFFFFFFFFu) table->total_regime_ticks[r]++;
entry->tick_count++;
entry->prev_window = current_window;
entry->prev_locality_q16 = current_locality_q16;
entry->has_prev = 1u;
/* --- UCB1: select best config for the next tick --- */
total = table->total_regime_ticks[r];
ln_total = cd_ln_approx_q16(total > 0u ? total : 1u);
best_c = 0u;
best_ucb = 0u;
for (ci = 0u; ci < table->num_configs; ci++) {
uint32_t nc = table->regime_ticks[r * table->num_configs + ci];
uint32_t ucb;
if (nc == 0u) {
ucb = tuning->score_max;
} else {
uint32_t ratio = ln_total / nc;
uint32_t sqrt_r = cd_isqrt32(ratio);
uint32_t bonus = (uint32_t)(((uint64_t)tuning->ucb_k * sqrt_r) / 256u);
uint32_t sc = table->regime_scores[r * table->num_configs + ci];
ucb = sc + bonus;
if (ucb > tuning->score_max) ucb = tuning->score_max;
}
if (ucb > best_ucb) { best_ucb = ucb; best_c = ci; }
}
table->current_config = (uint8_t)best_c;
return table->entries[best_c].config_bits;
}
void cd_config_table_force(CDConfigTable *table, uint8_t config_idx)
{
if (!table || !table->entries) return;
if ((uint32_t)config_idx >= table->num_configs) return;
table->current_config = config_idx;
}