-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathbasic_examples.py
More file actions
159 lines (121 loc) · 4.41 KB
/
Copy pathbasic_examples.py
File metadata and controls
159 lines (121 loc) · 4.41 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
"""
Basic examples of the Puiseux differential algebra pipeline.
Run with:
python examples/basic_examples.py
Requires:
- sympy
- DifferentialAlgebra (compiled C extension)
"""
from sympy import var
from DifferentialAlgebra import indexedbase
# Setup shared symbols
xi = var("xi")
x, y = indexedbase("x,y")
SEP = "=" * 70
def example_1_dy2_minus_x():
"""ẏ² - x = 0 → y = C + (2/3)x^{3/2}"""
print(f"\n{SEP}")
print("Example 1: ẏ² - x = 0")
print(f"Expected: y = z₀ + (2/3)x^{3/2}")
print(SEP)
from src.puiseux_pipeline import find_puiseux_solutions, filter_nontrivial
P = y[xi] ** 2 - x
sols = find_puiseux_solutions(P, b_max=10, verbose=True)
sols = filter_nontrivial(sols)
print(f"\nNon-trivial solutions found: {len(sols)}")
for s in sols:
print(f" (a,b) = ({s['a']}, {s['b']})")
print(f" FPS series : {s['series']}")
print(f" Puiseux : {s['puiseux_solution']}")
def example_2_dy3_minus_x2():
"""ẏ³ - x² = 0 → y = C + (3/5)x^{5/3}"""
print(f"\n{SEP}")
print("Example 2: ẏ³ - x² = 0")
print(f"Expected: y = z₀ + (3/5)x^{5/3}")
print(SEP)
from src.puiseux_pipeline import find_puiseux_solutions, filter_nontrivial
P = y[xi] ** 3 - x**2
sols = find_puiseux_solutions(P, b_max=15, verbose=True)
sols = filter_nontrivial(sols)
print(f"\nNon-trivial solutions found: {len(sols)}")
for s in sols:
print(f" (a,b) = ({s['a']}, {s['b']})")
print(f" Puiseux: {s['puiseux_solution']}")
def example_3_class_P_32_1():
"""8ẏ³ - 24ẏ² - 243y + 32 = 0 → y = 2x + 3x^{3/2}"""
print(f"\n{SEP}")
print("Example 3: Class P_{(3/2,1)}")
print("8ẏ³ - 24ẏ² - 243y + 32 = 0")
print(f"Expected: y = 2x + 3x^{3/2}")
print(SEP)
from src.puiseux_pipeline import find_puiseux_solutions, filter_nontrivial
P = 8 * y[xi] ** 3 - 24 * y[xi] ** 2 - 243 * y + 32
sols = find_puiseux_solutions(P, b_max=10, verbose=True)
sols = filter_nontrivial(sols)
print(f"\nNon-trivial solutions found: {len(sols)}")
for s in sols:
print(f" (a,b) = ({s['a']}, {s['b']})")
if s.get("puiseux_solution") is not None:
print(f" Puiseux: {s['puiseux_solution']}")
def example_4_fast_decision():
"""Demonstrate the FAST DECISION algorithm."""
print(f"\n{SEP}")
print("Example 4: Fast Decision Algorithm")
print(SEP)
from src.two_algorithms import fast_decision
test_cases = [
("ẏ² - x", y[xi] ** 2 - x),
("ẏ³ - x²", y[xi] ** 3 - x**2),
("8ẏ³ - 24ẏ² - 243y + 32",
8 * y[xi] ** 3 - 24 * y[xi] ** 2 - 243 * y + 32),
("ẏ³ - x", y[xi] ** 3 - x), # should be NO
]
for name, P in test_cases:
result, reason = fast_decision(P, verbose=False)
print(f"\n {name}:")
print(f" Decision: {result}")
print(f" Reason: {reason}")
def example_5_second_order():
"""ÿ² - y - x² - 1 = 0 (second-order AODE)"""
print(f"\n{SEP}")
print("Example 5: Second-order AODE")
print("ÿ² - y - x² - 1 = 0")
print(SEP)
from src.puiseux_pipeline import find_puiseux_solutions, filter_nontrivial
P = y[xi, xi] ** 2 - y - x**2 - 1
sols = find_puiseux_solutions(P, b_max=8, verbose=True)
sols = filter_nontrivial(sols)
print(f"\nNon-trivial solutions found: {len(sols)}")
for s in sols:
print(f" (a,b) = ({s['a']}, {s['b']})")
if s.get("puiseux_solution") is not None:
print(f" Puiseux: {s['puiseux_solution']}")
if __name__ == "__main__":
print("Puiseux Differential Algebra — Examples")
print("=" * 70)
print("Make sure DifferentialAlgebra is installed:")
print(" cd DifferentialAlgebra/bmi-install-sympy-c99/packaging")
print(" python3 setup.py build && pip install -e .")
print()
try:
example_1_dy2_minus_x()
except Exception as e:
print(f"Example 1 failed: {e}")
try:
example_2_dy3_minus_x2()
except Exception as e:
print(f"Example 2 failed: {e}")
try:
example_3_class_P_32_1()
except Exception as e:
print(f"Example 3 failed: {e}")
try:
example_4_fast_decision()
except Exception as e:
print(f"Example 4 failed: {e}")
try:
example_5_second_order()
except Exception as e:
print(f"Example 5 failed: {e}")
print(f"\n{SEP}")
print("Done.")