From 580776520f122ca9c20ac458206f771e36721a92 Mon Sep 17 00:00:00 2001
From: lhye200 <lhye20507@outlook.com>
Date: Wed, 23 Apr 2025 20:51:21 +0800
Subject: [PATCH] 20250423

---
 120-4.py | 25 +++++++++++++++++++++++++
 120-5.py | 55 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 80 insertions(+)
 create mode 100644 120-4.py
 create mode 100644 120-5.py

diff --git a/120-4.py b/120-4.py
new file mode 100644
index 0000000..c5f9e2c
--- /dev/null
+++ b/120-4.py
@@ -0,0 +1,25 @@
+def fx(x):
+    return 1/x
+
+
+def SplitTrapezoidal(a,b,err):
+    t1 = (b-a)*(fx(a)+fx(b))/2
+    k = 1
+    while True:
+        tmp = 0
+        for i in range(1, 2**(k-1)+1):
+            tmp += fx(a + (b-a)*(2*i-1)/(2**k))
+        t2 = t1/2+(b-a)*tmp/(2**k)
+        if abs(t2-t1) < err:
+            break
+        t1 = t2
+        k += 1
+    return t2,k
+
+
+if __name__ == "__main__":
+    a = 1
+    b = 3
+    err = 1e-2
+    result,k = SplitTrapezoidal(a, b, err)
+    print(f"Result: {result},k={k}")
diff --git a/120-5.py b/120-5.py
new file mode 100644
index 0000000..c4e37fe
--- /dev/null
+++ b/120-5.py
@@ -0,0 +1,55 @@
+import math
+
+def fx(x):
+    if x == 0:
+        x = 1e-10  # Avoid division by zero
+    return math.sin(x)/x
+
+
+
+def Romberg(a, b, err):
+    t00 = (b-a)*(fx(a)+fx(b))/2
+    t01 = 0
+    t10 = 0
+    t11 = 0
+    t20 = 0
+    t21 = 0
+    t30 = 0
+    t31 = 0
+    k = 1
+    while True:
+        tmp = 0
+        for i in range(1, 2**(k-1)+1):
+            tmp += fx(a + (b-a)*(2*i-1)/(2**k))
+        t01 = t00/2+(b-a)*tmp/(2**k)
+        
+        if k>1:
+            t11 = (4*t01-t00)/3
+            
+            if k>2:
+                t21 = (16*t11-t10)/15
+                
+                if k>3:
+                    t31 = (64*t21-t20)/63
+                    if abs(t31-t30) < err:
+                        break
+                    t30 = t31
+                else:
+                    t30 = (64*t21-t20)/63
+                t20 = t21
+            else:
+                t20 = (16*t11-t10)/15
+            t10 = t11
+        else:
+            t10 = (4*t01-t00)/3
+        t00 = t01
+        k += 1
+    return t31, k
+
+
+if __name__ == "__main__":
+    a = 0
+    b = 1
+    err = 0.5e-6
+    result, k = Romberg(a, b, err)
+    print(f"Result: {result}, k={k}")