import unittest
import sympy as s
from certificates import Interval,ResidualCertificate,BlockCandidates,polynomial_box
from models import Reaction,MassActionModel,LoadedReactor,EnvZOmpR,envz_box_floors,PrivateRelease
from invariant_region import InvariantEllipsoid
R=s.Rational

class ScientificChecks(unittest.TestCase):
    def test_candidates_are_not_acr_decisions(self):
        a,b,c=s.symbols('a b c');field=(b*(a-2),c*(a-3)+b)
        self.assertEqual(BlockCandidates((a,b,c),a).compute(field)['candidates'],['2','3'])
        u,v,t=s.symbols('u v t');self.assertEqual(BlockCandidates((u,v,t),t).compute((t*u-v,t*v-u))['candidates'],['1'])
    def test_signed_certificate_and_two_failure_types(self):
        model=LoadedReactor().model();a,b,c=model.variables;cert=ResidualCertificate(model.variables,0,R(11,5),b/20,(0,1,0));box=(Interval(0,10),Interval(1,5),Interval(0,10));eps=(0,R(1,1000),0)
        result=cert.evaluate(model.field,box,eps,R(1,20),(-a,0,0));self.assertEqual(result['absolute_error_bound'],'1/50');self.assertEqual(result['signed_load_bound'],'0')
        self.assertEqual(cert.evaluate(model.field,(box[0],Interval(0,5),box[2]),eps,R(1,20))['status'],'unresolved_multiplier_floor')
        self.assertEqual(ResidualCertificate(model.variables,0,2,b/20,(0,1,0)).evaluate(model.field,box,eps,R(1,20))['status'],'invalid_identity')
        with self.assertRaises(ValueError):cert.evaluate(model.field,box,(0,-1,0),R(1,20))
    def test_exact_operating_region_and_stability(self):
        e=LoadedReactor().equilibrium();self.assertEqual(e['state'],(R(11,5),R(17,10),R(17,10)));self.assertEqual(e['floor_load_limit'],R(29,1100));self.assertEqual(e['critical_load'],R(39,1100));self.assertTrue(e['locally_stable'])
        self.assertEqual(LoadedReactor(load=R(39,1100)).equilibrium()['status'],'no_positive_equilibrium');self.assertFalse(LoadedReactor(load=R(3,100)).equilibrium()['floor_satisfied'])
    def test_invariant_region_is_operating_point_specific(self):
        r=LoadedReactor();e=r.equilibrium();proof=InvariantEllipsoid(r.model(),e['state']).verify();self.assertEqual(proof['status'],'certified_forward_invariant');self.assertEqual(proof['multiplier_b_floor'],R(1573,1000))
        changed=LoadedReactor(load=R(1,100));self.assertEqual(InvariantEllipsoid(changed.model(),changed.equilibrium()['state']).verify()['status'],'not_applicable_to_this_operating_point')
        self.assertEqual(InvariantEllipsoid(r.model(),e['state'],10).verify()['status'],'not_certified_at_requested_level')
    def test_envz_full_field_pools_and_regularity(self):
        env=EnvZOmpR();model=env.model();eq=env.equilibrium(R(3,2),R(9,2));z=eq['state'];self.assertEqual(model.regularity(z)['status'],'regular');self.assertEqual(s.simplify(sum(z[i] for i in [0,1,2,5,6])),R(3,2));self.assertEqual(s.simplify(sum(z[i] for i in [3,4,5,6])),R(9,2));self.assertEqual(env.equilibrium(1,2)['status'],'no_positive_equilibrium')
        result=envz_box_floors((Interval('99/100','101/100'),)*9,Interval(1,2),Interval(4,5));self.assertEqual(result['x2_floor'],R(9121780899,77264239204))
    def test_private_release_and_transported_residual(self):
        old=MassActionModel(('a','b'),[Reaction((1,1),(0,3),1),Reaction((0,1),(1,0),1),Reaction((0,1),(0,0),1)]);new,W=PrivateRelease(0,(2,3)).apply(old);a,b,z1,z2=new.variables;sub={z1:a*b/2,z2:a*b/3};self.assertTrue(all(s.expand(new.field[i].subs(sub)-old.field[i])==0 for i in [0,1]));self.assertEqual(s.expand((a-1)*b-new.field[1]/2-R(3,2)*new.field[2]-new.field[3]),0)
        c=PrivateRelease(0,(1,1),('1/100','1/100')).steady_costs(1);self.assertEqual(c['product_yields'],[R(100,101),R(10000,10201)])
    def test_intervals_and_no_false_floor(self):
        x=s.Symbol('x');self.assertEqual(polynomial_box(x*x,(x,),(Interval(-1,2),)).strings(),['0','4'])
        # Exact polynomial simplification occurs before interval evaluation.
        self.assertEqual(polynomial_box(x-x,(x,),(Interval(-10,10),)).strings(),['0','0'])
        with self.assertRaises(ValueError):Interval(1,2)/Interval(-1,1)

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