SystemSimulationApp/tests/test_native_properties.py

"""Property reuse must preserve state identity, trial isolation and EOS physics."""
from pathlib import Path
import subprocess
import tempfile
import unittest

from app.simulation.native_codegen.build import toolchain
from tests.native_kernel_source import native_kernel_source

ROOT = Path(__file__).resolve().parents[1]


class NativePropertyTests(unittest.TestCase):
    def test_shared_properties_and_invalidation(self):
        try:
            compiler = toolchain()[0]
        except (OSError, RuntimeError, subprocess.SubprocessError) as exc:
            self.skipTest(f"Native toolchain unavailable: {exc}")
        source = native_kernel_source()
        for signature, counter in (
            ('static double temperature_ph(double p,double h) {', 'ph_calls'),
            ('static double z_factor(double p,double T) {', 'z_calls'),
            ('double native_viscosity(const NativeMedium *m, double T, int diagnostic) {', 'mu_calls'),
        ):
            self.assertEqual(source.count(signature), 1)
            source = source.replace(signature, signature + f'++{counter};')
        source = 'static int ph_calls,z_calls,mu_calls;\n' + source
        harness = r'''
#include <stdio.h>
#define CHECK(x) do { if(!(x)){fprintf(stderr,"line %d\n",__LINE__);return 1;} } while(0)
static int close_to(double a,double b,double atol,double rtol) {
    return isfinite(a) && isfinite(b) && fabs(a-b)<=atol+rtol*fmax(fabs(a),fabs(b));
}
int main(void) {
    NativeMedium m=helium_medium;
    NativePropertyState states[32],other_states[2];NativePropertyCache cache,other;
    native_properties_init(&cache,states,32);native_properties_init(&other,other_states,2);
    double y[2];NativeGas gas;
    CHECK(native_medium_init(&m,1.5e7,293.15,.01,0,y));
    int before=z_calls;
    CHECK(native_medium_gas_context(&cache,&m,y[0],y[1],.01,&gas));
    CHECK(z_calls==before); /* h=u+p/rho does not solve the cubic again. */
    CHECK(gas.h==gas.u+gas.p/gas.rho);
    before=ph_calls;
    CHECK(native_temperature_ph_context(&cache,&m,gas.p,gas.h)==gas.T && ph_calls==before);
    NativeMedium copy=m;
    CHECK(native_temperature_ph_context(&cache,&copy,gas.p,gas.h)==gas.T && ph_calls==before);
    double q,cm,v;before=z_calls;
    CHECK(native_medium_orifice_context(&cache,&m,gas.p,.6*gas.p,gas.h,gas.h,.001,.5,&q,&cm,&v));
    CHECK(z_calls==before+1); /* Only the different downstream state needs rho. */
    CHECK(native_medium_orifice_context(&cache,&m,gas.p,.6*gas.p,gas.h,gas.h,.001,.5,&q,&cm,&v));
    CHECK(z_calls==before+1);
    before=mu_calls;
    q=native_pipe_flow_context(&cache,&m,gas.p,.6*gas.p,gas.T,.01,1,1e-5,1);
    CHECK(isfinite(q) && mu_calls==before+1);
    double diag[4];native_pipe_diagnostics_context(&cache,&m,q,gas.p,gas.T,.01,1,1e-5,1,diag);
    CHECK(mu_calls==before+1);

    /* Changed pressure (even one ULP), mixed enthalpy and every medium field
       must force a fresh PH evaluation. Equal h alone is not a state key. */
    double pressures[]={nextafter(gas.p,INFINITY),.2*gas.p,gas.p};
    for(int i=0;i<3;i++) {
        double h=gas.h+(i==2?1000:0),expected=native_temperature_ph(&m,pressures[i],h);
        before=ph_calls;
        CHECK(native_temperature_ph_context(&cache,&m,pressures[i],h)==expected);
        CHECK(ph_calls==before+1);
        CHECK(native_temperature_ph_context(&cache,&m,pressures[i],h)==expected && ph_calls==before+1);
    }
    CHECK(fabs(native_temperature_ph(&m,.2*gas.p,gas.h)-gas.T)>1e-3);
    double positive_h=gas.h+2e6;
    (void)native_temperature_ph_context(&cache,&m,gas.p,positive_h);
    for(int i=0;i<8;i++) {
        copy=m;
        double *fields[]={&copy.R,&copy.cp,&copy.Tref,&copy.slope,&copy.mu,&copy.muT,&copy.S};
        if(i==7)copy.real_helium=0;else *fields[i]+=fmax(fabs(*fields[i])*.001,.00001);
        double expected=native_temperature_ph(&copy,gas.p,positive_h);
        size_t count=cache.count;
        CHECK(native_temperature_ph_context(&cache,&copy,gas.p,positive_h)==expected);
        CHECK(cache.count>count);
    }
    /* Separate contexts, a new trial, and exhausted capacity remain correct. */
    before=ph_calls;
    CHECK(isfinite(native_temperature_ph_context(&other,&m,gas.p,gas.h)) && ph_calls==before+1);
    native_properties_init(&cache,states,1);
    CHECK(native_medium_gas_context(&cache,&m,y[0],y[1],.01,&gas));
    for(int i=1;i<30;i++) {
        double p=gas.p*(1+i*.01),h=gas.h+i;
        CHECK(native_temperature_ph_context(&cache,&m,p,h)==native_temperature_ph(&m,p,h));
        CHECK(cache.count==1);
    }
    native_properties_init(&cache,states,32);before=ph_calls;
    CHECK(isfinite(native_temperature_ph_context(&cache,&m,gas.p,gas.h)) && ph_calls==before+1);
    size_t count=cache.count;
    (void)native_temperature_ph_context(&cache,&m,NAN,gas.h);
    CHECK(cache.count==count);

    /* Thermodynamic round trips use an independent PT enthalpy formula.
       Ideal gases include temperature-dependent heat capacity. */
    for(int ideal=0;ideal<2;ideal++) {
        NativeMedium fluid=ideal?(NativeMedium){0,287,1005,300,.2,1.8e-5,300,110.4}:m;
        for(double T=100;T<=1000;T+=100)for(double p=1e4;p<=1e8;p*=10) {
            native_properties_init(&cache,states,32);
            CHECK(native_medium_init(&fluid,p,T,.01,0,y));
            CHECK(native_medium_gas_context(&cache,&fluid,y[0],y[1],.01,&gas));
            double ref_h=ideal?fluid.cp*T+.5*fluid.slope*pow(T-fluid.Tref,2):h_ideal(T)+h_departure(p,T);
            CHECK(close_to(gas.h,ref_h,2e-5,1e-9));
            CHECK(close_to(gas.p,p,.001,1e-9));
            CHECK(close_to(gas.T,T,1e-7,1e-9));
            CHECK(close_to(native_temperature_ph(&fluid,gas.p,gas.h),gas.T,1e-7,1e-9));
            double oldq,oldcm,oldv;
            CHECK(native_medium_orifice(&fluid,gas.p,.8*gas.p,gas.h,gas.h,.001,.5,&oldq,&oldcm,&oldv));
            CHECK(native_medium_orifice_context(&cache,&fluid,gas.p,.8*gas.p,gas.h,gas.h,.001,.5,&q,&cm,&v));
            CHECK(close_to(q,oldq,1e-10,2e-9) && close_to(v,oldv,1e-7,2e-9));
        }
    }
    /* Subcritical helium retains the pre-existing vapor-root selection. */
    native_properties_init(&cache,states,32);
    CHECK(native_gas_init(1e4,3,.01,y));
    CHECK(native_medium_gas_context(&cache,&m,y[0],y[1],.01,&gas));
    CHECK(gas.T<TC && cache.count==0);
    CHECK(gas.h==h_ideal(gas.T)+h_departure(gas.p,gas.T));
    return 0;
}
'''
        with tempfile.TemporaryDirectory(prefix='native-properties-') as tmp:
            folder=Path(tmp);c_file=folder/'properties.c';exe=folder/'properties.exe'
            c_file.write_text(source+harness)
            compiled=subprocess.run([compiler,'-std=c11','-O3','-Wall','-Wextra','-Werror',
                '-ffp-contract=off','-fno-fast-math','-static-libgcc','-I',str(ROOT/'native/include'),
                str(c_file),'-lm','-o',str(exe)],capture_output=True,text=True,timeout=60)
            self.assertEqual(compiled.returncode,0,compiled.stderr)
            run=subprocess.run([str(exe)],capture_output=True,text=True,timeout=30)
            self.assertEqual(run.returncode,0,run.stderr)


if __name__ == '__main__':
    unittest.main()