Light Transport

#include "lantern.hpp"
#include <cmath>
#include <iostream>
#include <limits>
#include <numbers>
#include <stdexcept>

namespace {
std::size_t checks{};
void require(bool condition,const char* message) { ++checks; if (!condition) throw std::runtime_error{message}; }
bool near(double a,double b,double error=1e-10) { return std::abs(a-b)<error; }
void sampling() {
    lantern::Rng known{0}; require(known.next()==0xe220a8397b1dcdafULL,"RNG known seed/word contract");
    lantern::Rng a{47},b{47};
    bool repeat=true,range=true;
    for (unsigned i=0;i<10000;++i) { const auto x=a.uniform(),y=b.uniform(); repeat &= x==y; range &= x>=0 && x<1; }
    require(repeat && range,"fixed seed repeats and uniform values stay half-open");
    const auto pole=lantern::cosine_hemisphere(0,0);
    require(pole && near(pole->direction.z,1) && near(pole->pdf,1/std::numbers::pi),"cosine sampler pole and steradian PDF");
    require(!lantern::cosine_hemisphere(1,0),"exclusive sampler endpoint is rejected");
    require(!lantern::cosine_hemisphere(0,std::numeric_limits<double>::quiet_NaN()),"NaN sampler input rejected");
    lantern::Rng rng{8}; double mean_z=0; bool valid=true;
    for (unsigned i=0;i<8192;++i) {
        const double u=rng.uniform(), v=rng.uniform(); // Fix sample dimension order.
        const auto s=lantern::cosine_hemisphere(u,v);
        if (!s) { valid=false; break; }
        const auto d=s->direction;
        valid &= near(d.x*d.x+d.y*d.y+d.z*d.z,1) && d.z>0 && near(s->pdf,d.z/std::numbers::pi);
        mean_z+=d.z;
    }
    require(valid && std::abs(mean_z/8192-2.0/3)<0.02,"sample directions and a known cosine-distribution moment");
    const auto estimate=lantern::cosine_integral(65536,42);
    require(estimate && near(*estimate,std::numbers::pi,0.03),"Monte Carlo estimate agrees with analytic hemisphere cosine integral");
    require(!lantern::cosine_integral(0,1) && !lantern::cosine_integral(lantern::max_samples+1,1),"sample work is bounded");
    const auto weight=lantern::power_weight(2,1), reverse=lantern::power_weight(1,2);
    require(weight && reverse && near(*weight,0.8) && near(*weight+*reverse,1),"power weights partition the two-strategy contribution");
    const auto tiny=lantern::power_weight(1e-300,1e-300);
    require(tiny && *tiny==0.5,"tiny PDFs do not underflow the MIS denominator");
    const auto subnormal=lantern::power_weight(std::numeric_limits<double>::denorm_min(),
                                              std::numeric_limits<double>::denorm_min());
    require(subnormal && *subnormal==0.5,"equal subnormal PDFs preserve the power weight");
    const auto absent_a=lantern::power_weight(0,1), absent_b=lantern::power_weight(1,0);
    require(absent_a && absent_b && *absent_a==0 && *absent_b==1,
            "zero-support strategies receive zero weight");
    require(!lantern::power_weight(0,0) && !lantern::power_weight(-1,1),"invalid PDF pairs rejected");
    lantern::Moments m;
    require(m.add(1).has_value() && m.add(2).has_value() && m.add(3).has_value(),"finite samples accumulate");
    require(m.count()==3 && near(m.mean(),2) && near(m.sample_variance(),1),"Welford mean and sample variance known values");
    require(!m.add(std::numeric_limits<double>::infinity()) && m.count()==3,"invalid sample leaves prior statistics intact");
}
void boxes() {
    const lantern::Aabb box{{-1,-1,-1},{1,1,1}};
    const auto central=lantern::box_hit({{0,0,4},{0,0,-1}},box);
    require(central && *central,"parallel-axis central box hit");
    const auto outside=lantern::box_hit({{2,0,4},{-0.0,0,-1}},box);
    require(outside && !*outside,"signed-zero direction outside slab misses without division");
    const auto tangent=lantern::box_hit({{1,0,4},{0,0,-1}},box);
    require(tangent && *tangent,"tangent slab boundary is accepted");
    const auto behind=lantern::box_hit({{0,0,4},{0,0,1}},box);
    require(behind && !*behind,"box behind ray is rejected");
    require(!lantern::box_hit({{0,0,4},{0,0,-2}},box),"nonunit ray contract checked");
    require(!lantern::box_hit({{0,0,4},{0,0,-1}},{{1,1,1},{-1,-1,-1}}),"inverted bounds rejected");
}
void acceleration() {
    auto empty=lantern::SceneIndex::build({});
    require(empty && empty->node_count()==0 && empty->closest({{0,0,4},{0,0,-1}})->has_value()==false,"empty index queries safely miss");
    std::vector<lantern::Sphere> input;
    for (unsigned i=0;i<64;++i) input.push_back({i,{static_cast<double>(i%8)-3.5,static_cast<double>(i/8)-3.5,0},0.35});
    auto built=lantern::SceneIndex::build(input);
    require(built && built->node_count()<=2*input.size()-1,"binary index uses bounded owned nodes");
    std::vector<lantern::Sphere>{}.swap(input); // Release caller storage, not just its size.
    lantern::Rng rng{100}; bool agreement=true;
    for (unsigned i=0;i<2048;++i) {
        const lantern::Ray ray{{rng.uniform()*10-5,rng.uniform()*10-5,5},{0,0,-1}};
        const auto tree=built->closest(ray),linear=built->linear_closest(ray);
        if (!tree || !linear || tree->has_value()!=linear->has_value() ||
            (*tree && ((*tree)->id!=(*linear)->id || !near((*tree)->distance,(*linear)->distance)))) { agreement=false; break; }
    }
    require(agreement,"BVH matches linear oracle over deterministic hit/miss rays after input owner cleared");
    std::vector<lantern::Sphere> tied{{9,{0,0,0},1},{2,{0,0,0},1},{4,{0,0,0},0.5}};
    const auto tie=lantern::SceneIndex::build(tied);
    const auto hit=tie->closest({{0,0,4},{0,0,-1}});
    require(hit && *hit && (*hit)->id==2 && near((*hit)->distance,3),"tied centroid build terminates and equal hit distances use smallest id");
    const auto inside=tie->closest({{0,0,0},{0,0,1}});
    require(inside && *inside && (*inside)->id==4 && near((*inside)->distance,0.5),"inside-sphere positive exit root preserved by BVH");
    const auto oblique=tie->closest({{3,4,0},{-0.6,-0.8,0}});
    require(oblique && *oblique && (*oblique)->id==2 && near((*oblique)->distance,4),
            "oblique central ray has independently known distance and tie id");
    std::vector<lantern::Sphere> maximum;
    for (std::size_t i=0; i<lantern::max_spheres; ++i)
        maximum.push_back({static_cast<std::uint32_t>(i), {0,0,0}, 1});
    const auto maximum_index=lantern::SceneIndex::build(maximum);
    require(maximum_index && maximum_index->node_count()<=2*lantern::max_spheres-1,
            "exact primitive budget builds with coincident centroids");
    const auto maximum_hit=maximum_index->closest({{0,0,4},{0,0,-1}});
    require(maximum_hit && *maximum_hit && (*maximum_hit)->id==0 && near((*maximum_hit)->distance,3),
            "exact-budget traversal retains the smallest tied id");
    const std::vector<lantern::Sphere> tangent_sphere{{1,{0,0,0},1}};
    const auto grazing=lantern::SceneIndex::build(tangent_sphere);
    for (const double x : {std::nextafter(1.0,2.0),1.0,std::nextafter(1.0,0.0)}) {
        const lantern::Ray ray{{x,0,1000},{0,0,-1}};
        const auto tree=grazing->closest(ray),linear=grazing->linear_closest(ray);
        const bool should_hit=x<=1;
        require(tree && linear && tree->has_value()==should_hit && linear->has_value()==should_hit,
                "far-origin tangent and adjacent representable rays agree with the geometric oracle");
    }
    tied[1].id=9; require(!lantern::SceneIndex::build(tied),"duplicate ids rejected");
    tied.resize(lantern::max_spheres+1); require(!lantern::SceneIndex::build(tied),"primitive budget rejected before copying");
    tied={{1,{0,0,0},std::numeric_limits<double>::quiet_NaN()}};
    require(!lantern::SceneIndex::build(tied),"non-finite primitive radius rejected");
    require(!built->closest({{0,0,0},{0,0,0}}),"invalid query rejected without traversal");
}
}
int main() {
    try { sampling(); boxes(); acceleration(); std::cout << checks << " Lantern CPU checks passed.\n"; return 0; }
    catch (const std::exception& error) { std::cerr << error.what() << '\n'; return 1; }
}