From 591495159ef29a3edac10aaf7e3d7ea28a4a3013 Mon Sep 17 00:00:00 2001 From: Deniz Donmez Date: Mon, 15 Sep 2025 10:40:52 +0300 Subject: [PATCH] new examples --- examples/demo.py | 91 +++++++++++++++++++++++++++++++++--------------- 1 file changed, 62 insertions(+), 29 deletions(-) diff --git a/examples/demo.py b/examples/demo.py index 5a19eff..a57f143 100644 --- a/examples/demo.py +++ b/examples/demo.py @@ -2,27 +2,6 @@ import sys sys.path.append('../') from numethods import * -A = Matrix([[2, -1], [1, 2], [1, 1]]) -b = Vector([1, 2, 3]) - -# Factorization -qr = QRHouseholder(A) -Q, R = qr.Q, qr.R -print("Q =", Q) -print("R =", R) - -qrm = QRModifiedGramSchmidt(A) -Qm, Rm = qrm.Q, qrm.R -print("Qm =", Qm) -print("Rm =", Rm) -print("Q^T Q =", Q.T @ Q) -print("Qm^T Qm =", Qm.T @ Qm) -print("A=Qm Rm =", Qm @ Rm) - -# Solve Ax = b (least squares, since A is tall) -x_ls = LeastSquaresSolver(A, b).solve() -print("Least squares solution:", x_ls) - def demo_linear_solvers(): A = Matrix([[4, -1, 0], [-1, 4, -1], [0, -1, 3]]) @@ -34,16 +13,56 @@ def demo_linear_solvers(): print("Jacobi:", Jacobi(A, b, tol=1e-12).solve()) print("Gauss-Seidel:", GaussSeidel(A, b, tol=1e-12).solve()) +def demo_qr(): + A = Matrix([[2, -1], [1, 2], [1, 1]]) + + b = Vector([1, 2, 3]) + + # Factorization + qr = QRHouseholder(A) + Q, R = qr.Q, qr.R + print("Q =", Q) + print("R =", R) + + qrm = QRModifiedGramSchmidt(A) + Qm, Rm = qrm.Q, qrm.R + print("Qm =", Qm) + print("Rm =", Rm) + print("Q^T Q =", Q.T @ Q) + print("Qm^T Qm =", Qm.T @ Qm) + print("A=Qm Rm =", Qm @ Rm) + print("A=Q R =", Q @ R) + + # Solve Ax = b (least squares, since A is tall) + x_ls = LeastSquaresSolver(A, b).solve() + print("Least squares solution:", x_ls) + + def demo_roots(): - f = lambda x: x**3 - x - 2 - df = lambda x: 3 * x**2 - 1 - print("Bisection:", Bisection(f, 1, 2).solve()) - print("Secant:", Secant(f, 1.0, 2.0).solve()) - print("Newton root:", NewtonRoot(f, df, 1.5).solve()) - # A simple contraction for demonstration; not general-purpose - g = lambda x: (x + 2 / x**2) / 2 - print("Fixed point (demo):", FixedPoint(g, 1.5).solve()) + f = lambda x: x**2 - 2 + df = lambda x: 2 * x + + # Newton + steps = NewtonRoot(f, df, x0=1.0).trace() + print("Newton Method Trace (x^2 - 2):") + print_trace(steps) + + # Secant + steps = Secant(f, 0, 2).trace() + print("\nSecant Method Trace (x^2 - 2):") + print_trace(steps) + + # Bisection + steps = Bisection(f, 0, 2).trace() + print("\nBisection Method Trace (x^2 - 2):") + print_trace(steps) + + # Fixed-point: solve + g = lambda x: 0.5 * (x + 2 / x) + steps = FixedPoint(g, 1.0).trace() + print("\nFixed-Point Iteration Trace (x^2 - 2):") + print_trace(steps) def demo_interpolation(): @@ -55,8 +74,22 @@ def demo_interpolation(): print("Newton interpolation at", t, "=", newt.evaluate(t)) print("Lagrange interpolation at", t, "=", lagr.evaluate(t)) +def demo_differentiation(): + f = lambda x: x**3 # f'(x) = 3x^2, f''(x) = 6x + x0 = 2.0 + + print("Forward :", ForwardDiff(f, x0)) + print("Backward :", BackwardDiff(f, x0)) + print("Central :", CentralDiff(f, x0)) + print("4th order:", CentralDiff4th(f, x0)) + print("Richardson:", RichardsonExtrap(f, x0)) + print("Second derivative:", SecondDerivative(f, x0)) + + if __name__ == "__main__": demo_linear_solvers() demo_roots() demo_interpolation() + demo_qr() + demo_differentiation()