Language: C++

// Template
#include "bits/stdc++.h"
#define rep_override(x, y, z, name, ...) name
#define rep2(i, n) for (int i = 0; i < (int)(n); ++i)
#define rep3(i, l, r) for (int i = (int)(l); i < (int)(r); ++i)
#define rep(...) rep_override(__VA_ARGS__, rep3, rep2)(__VA_ARGS__)
#define per(i, n) for (ll i = (ll)(n) - 1; i >= 0; --i)
#define all(x) (x).begin(), (x).end()
using namespace std;
using ll = long long;
constexpr int inf = 1001001001;
constexpr ll INF = 3003003003003003003LL;
template <typename T>
inline bool chmin(T &x, const T &y) {
    if (x > y) {
        x = y;
        return true;
    }
    return false;
}
template <typename T> 
inline bool chmax(T &x, const T &y) {
    if (x < y) {
        x = y;
        return true;
    }
    return false;
}
struct IOSET {
    IOSET() {
        cin.tie(nullptr);
        ios::sync_with_stdio(false);
        cout << fixed << setprecision(10);
    }
} ioset;
struct TIMER {
    std::chrono::high_resolution_clock::time_point st;
    TIMER() {reset();}
    void reset() {st = std::chrono::high_resolution_clock::now();}
    std::chrono::milliseconds::rep get() {
        auto ed = std::chrono::high_resolution_clock::now();
        return std::chrono::duration_cast<std::chrono::milliseconds>(ed - st).count();
    }
} timer;
struct RANDOM {
    std::mt19937 mt;
    RANDOM() {
        std::random_device rd;
        mt = std::mt19937(rd());
    }
    int operator()(int a, int b) {
        std::uniform_int_distribution<int> dist(a, b - 1);
        return dist(mt);
    }
    int operator()(int a) {
        return (*this)(0, a);
    }
    bool prob(double x) {
        std::uniform_int_distribution<double> dist(0, 1);
        return dist(mt) < x;
    }
} rnd;

// Main
int main() {
    int d;
    cin >> d;
    vector<int> c(26);
    rep(i, 26) cin >> c[i];
    vector<vector<int>> s(d, vector<int>(26));
    rep(i, d) rep(j, 26) cin >> s[i][j];

    vector<int> t(d);
    vector<vector<int>> ds(26);
    int score = 0;
    auto rescore = [&](int x, int y) {
        if (y == t[x]) return;
        int mn = -inf, mx = inf;
        for (int i : ds[t[x]]) {
            if (i < x) chmax(mn, i);
            if (i > x) chmin(mx, i);
        }
        score -= s[x][t[x]];
        score += (x - mn) * (x - mn - 1) / 2 * c[t[x]];
        score += (mx - x) * (mx - x - 1) / 2 * c[t[x]];
        score -= (mx - mn) * (mx - mn - 1) / 2 * c[t[x]];
        mn = -inf, mx = inf;
        for (int i : ds[y]) {
            if (i < x) chmax(mn, i);
            if (i > x) chmin(mx, i);
        }
        score += s[x][y];
        score += (mx - mn) * (mx - mn - 1) / 2 * c[y];
        score -= (x - mn) * (x - mn - 1) / 2 * c[y];
        score -= (mx - x) * (mx - x - 1) / 2 * c[y];
        vector<int> new_;
        for (int i : ds[t[x]]) if (i != x) new_.push_back(i);
        ds[t[x]] = new_;
        new_.clear();
        rep(i, ds[y].size()) {
            if (ds[y][i] > x && ds[y][i - 1] < x) new_.push_back(x);
            new_.push_back(ds[y][i]);
        }
        ds[y] = new_;
        t[x] = y;
    };

    {
    vector<int> last(26, -1);
    rep(i, d) {
        int mx = -inf, idx = -1;
        rep(j, 26) {
            if (chmax(mx, (i - last[j]) * c[j] * 13 + s[i][j])) idx = j;
        }
        t[i] = idx;
        last[t[i]] = i;
    }
    last = vector<int>(26, -1);
    rep(i, d) {
        score += s[i][t[i]];
        last[t[i]] = i;
        rep(j, 26) score -= (i - last[j]) * c[j];
    }
    rep(i, 26) ds[i].push_back(-1);
    rep(i, d) ds[t[i]].push_back(i);
    rep(i, 26) ds[i].push_back(d);
    }

    double T0 = 1000.0, T1 = 100.0, T = T0;
    int cnt = 0;
    while (timer.get() < 1950) {
        ++cnt;
        if (cnt % 100 == 0) {
            double x = timer.get() / 1950.0;
            T = pow(T0, x) * pow(T1, 1 - x);
        }
        if (rnd.prob(0.5)) {
            int x = rnd(d), y = rnd(26), old = score, memo = t[x];
            rescore(x, y);
            if (score < old && !rnd.prob(exp((old - score) / T))) rescore(x, memo);
        } else {
            int x = rnd(d), y = rnd(min(x + 20, d)), memox = t[x], memoy = t[y], old = score;
            rescore(x, memoy);
            rescore(y, memox);
            if (score < old && !rnd.prob((old - score) / T)) {
                rescore(x, memox);
                rescore(y, memoy);
            }
        }
    }

    rep(i, d) cout << t[i] + 1 << "\n";
    cout << score << "\n";
}

Compiler output:
In file included from /usr/include/c++/14/bits/stl_algo.h:65,
                 from /usr/include/c++/14/algorithm:61,
                 from /usr/include/x86_64-linux-gnu/c++/14/bits/stdc++.h:51,
                 from a.cc:2:
/usr/include/c++/14/bits/uniform_int_dist.h: In instantiation of 'class std::uniform_int_distribution<double>':
a.cc:59:51:   required from here
   59 |         std::uniform_int_distribution<double> dist(0, 1);
      |                                                   ^
/usr/include/c++/14/bits/uniform_int_dist.h:90:49: error: static assertion failed: template argument must be an integral type
   90 |       static_assert(std::is_integral<_IntType>::value,
      |                                                 ^~~~~
/usr/include/c++/14/bits/uniform_int_dist.h:90:49: note: 'std::integral_constant<bool, false>::value' evaluates to false