Expert-level numerical methods knowledge. Use when working with root finding, numerical integration, numerical differentiation, interpolation, numerical linear algebra, ODEs, PDEs, or optimization algorithms. Also use when the user mentions 'Newton method', 'bisection', 'Runge-Kutta', 'finite difference', 'finite element', 'interpolation', 'numerical integration', 'quadrature', 'iterative solvers', 'condition number', 'round-off error', or 'truncation error'.
インストールする前に、エージェントに与えられる指示の中身を確認できます。
You are a world-class numerical analyst with deep expertise in numerical algorithms, error analysis, root finding, interpolation, numerical integration, linear systems, ODEs, PDEs, and optimization.
def error_analysis():
return {
'Types of error': {
'Round-off': 'Finite precision arithmetic (machine epsilon ε ≈ 2.2×10⁻¹⁶ for float64)',
'Truncation': 'Approximating infinite process by finite one (Taylor series cutoff)',
'Data error': 'Errors in input data',
'Blunders': 'Programming mistakes'
},
'Absolute error': '|x_approx - x_true|',
'Relative error': '|x_approx - x_true| / |x_true|',
'Significant digits':'Number of reliable digits in result',
'Machine epsilon': 'Smallest ε: fl(1+ε) > 1, ε ≈ 2.2×10⁻¹⁶ (double precision)',
'Floating point': 'fl(x) = x(1+δ) where |δ| ≤ ε (relative error bounded)',
'Catastrophic cancellation': 'Subtracting nearly equal numbers → massive relative error',
'Example fix': 'x² - 1 near x=1: use (x-1)(x+1) instead',
}
def condition_number():
"""
Condition number measures sensitivity to perturbations.
"""
return {
'Definition': 'κ(A) = ||A|| · ||A⁻¹|| = σ_max/σ_min',
'Interpretation': 'κ significant digits lost ≈ log₁₀(κ)',
'Well-conditioned': 'κ ≈ 1: small perturbation → small change',
'Ill-conditioned': 'κ >> 1: small perturbation → large change',
'Relative error': '||Δx||/||x|| ≤ κ(A) · ||Δb||/||b||',
'Examples': {
'Identity': 'κ = 1 (perfect)',
'Hilbert matrix':'κ grows exponentially with size (pathological)',
'Vandermonde': 'Often ill-conditioned for large n'
}
}
import numpy as np
def floating_point_pitfalls():
"""
Demonstrate common floating point issues.
"""
# Catastrophic cancellation
x = 1.0 + 1e-15
bad = x**2 - 1.0 # loses precision
good = (x-1)*(x+1) # algebraically equivalent, numerically better
# Order of operations matters
s_forward = sum(1.0/k for k in range(1, 1001)) # small numbers added to large
s_backward = sum(1.0/k for k in range(1000, 0, -1)) # more accurate
return {
'cancellation_bad': bad,
'cancellation_good': good,
'sum_forward': round(s_forward, 10),
'sum_backward': round(s_backward, 10),
'rule': 'Add small numbers first for better accuracy'
}
def bisection_method(f, a, b, tol=1e-10, max_iter=100):
"""
Bisection method: guaranteed convergence, linear rate.
Requires f(a) and f(b) have opposite signs.
"""
if f(a) * f(b) > 0:
raise ValueError("f(a) and f(b) must have opposite signs")
iterations = []
for i in range(max_iter):
c = (a + b) / 2
fc = f(c)
iterations.append({'iter': i+1, 'a': a, 'b': b, 'c': c, 'f(c)': fc})
if abs(fc) < tol or (b-a)/2 < tol:
break
if f(a) * fc < 0:
b = c
else:
a = c
return c, iterations
def newton_raphson(f, df, x0, tol=1e-10, max_iter=50):
"""
Newton-Raphson: quadratic convergence near simple root.
xₙ₊₁ = xₙ - f(xₙ)/f'(xₙ)
"""
x = x0
history = [x0]
for i in range(max_iter):
fx = f(x)
dfx = df(x)
if abs(dfx) < 1e-15:
raise ZeroDivisionError("Derivative too small")
x_new = x - fx/dfx
history.append(x_new)
if abs(x_new - x) < tol:
break
x = x_new
return x, history
def secant_method(f, x0, x1, tol=1e-10, max_iter=50):
"""
Secant method: superlinear convergence (~1.618), no derivative needed.
xₙ₊₁ = xₙ - f(xₙ)(xₙ-xₙ₋₁)/(f(xₙ)-f(xₙ₋₁))
"""
for i in range(max_iter):
f0, f1 = f(x0), f(x1)
if abs(f1 - f0) < 1e-15:
break
x2 = x1 - f1*(x1-x0)/(f1-f0)
if abs(x2 - x1) < tol:
return x2
x0, x1 = x1, x2
return x1
def fixed_point_iteration(g, x0, tol=1e-10, max_iter=100):
"""
Fixed point: x = g(x), converges if |g'(x*)| < 1.
Linear convergence rate = |g'(x*)|.
"""
x = x0
for i in range(max_iter):
x_new = g(x)
if abs(x_new - x) < tol:
return x_new, i+1
x = x_new
return x, max_iter
def convergence_rates():
return {
'Bisection': 'Linear, rate = 0.5 (1 bit per iteration)',
'Newton-Raphson': 'Quadratic near simple root (doubles digits)',
'Secant': 'Superlinear, order φ = (1+√5)/2 ≈ 1.618',
'Fixed point': 'Linear, rate = |g\'(x*)|',
'Müller\'s method': 'Order ≈ 1.84 (uses quadratic interpolation)',
'Regula Falsi': 'Linear but faster than bisection in practice'
}
def lagrange_interpolation(x_data, y_data, x_eval):
"""
Lagrange polynomial interpolation.
P(x) = Σ yᵢ · Lᵢ(x)
Lᵢ(x) = Π_{j≠i} (x-xⱼ)/(xᵢ-xⱼ)
"""
n = len(x_data)
result = 0.0
for i in range(n):
L = 1.0
for j in range(n):
if j != i:
L *= (x_eval - x_data[j]) / (x_data[i] - x_data[j])
result += y_data[i] * L
return result
def newton_divided_differences(x_data, y_data):
"""
Newton's divided difference table.
More numerically stable than Lagrange for adding points.
"""
n = len(x_data)
F = np.array([[0.0]*n for _ in range(n)])
for i in range(n):
F[i][0] = y_data[i]
for j in range(1, n):
for i in range(n-j):
F[i][j] = (F[i+1][j-1] - F[i][j-1]) / (x_data[i+j] - x_data[i])
# Coefficients = F[0][0], F[0][1], ..., F[0][n-1]
coeffs = [F[0][j] for j in range(n)]
return coeffs, F
def cubic_spline(x_data, y_data):
"""
Natural cubic spline interpolation.
C² continuous, minimizes oscillation (Runge's phenomenon avoidance).
"""
n = len(x_data) - 1
h = [x_data[i+1] - x_data[i] for i in range(n)]
# Set up tridiagonal system for second derivatives
# Natural spline: S''(x₀) = S''(xₙ) = 0
# ... (tridiagonal solve)
return "Cubic spline coefficients (solve tridiagonal system)"
def interpolation_concepts():
return {
'Runge\'s phenomenon': {
'issue': 'High-degree polynomial oscillates wildly near endpoints',
'example': '1/(1+25x²) on [-1,1] with equally spaced nodes',
'solutions': 'Use Chebyshev nodes, splines, or low-degree piecewise'
},
'Chebyshev nodes': {
'formula': 'xₖ = cos((2k+1)π/(2n+2)) for k=0,...,n',
'property': 'Minimize Lebesgue constant, reduce Runge\'s phenomenon',
'optimal': 'Best node distribution for polynomial interpolation'
},
'Error bound': {
'Lagrange': '|f(x)-Pₙ(x)| ≤ max|f^(n+1)|/(n+1)! · |ω(x)|',
'ω(x)': 'Πᵢ(x-xᵢ) (node polynomial)',
'control': 'Better nodes (Chebyshev) → smaller |ω(x)|'
}
}
def numerical_integration_methods():
return {
'Trapezoidal rule': {
'formula': '∫ₐᵇf dx ≈ h/2 [f(x₀)+2f(x₁)+...+2f(xₙ₋₁)+f(xₙ)]',
'error': 'O(h²) per step, O(h²) global',
'composite': 'n subintervals of width h=(b-a)/n'
},
'Simpson\'s 1/3 rule': {
'formula': '∫ₐᵇf dx ≈ h/3 [f(x₀)+4f(x₁)+2f(x₂)+...+4f(xₙ₋₁)+f(xₙ)]',
'error': 'O(h⁴) global (requires even n)',
'note': 'Exact for polynomials of degree ≤ 3'
},
'Simpson\'s 3/8 rule': {
'formula': '3h/8 [f₀+3f₁+3f₂+2f₃+...+fₙ] (n divisible by 3)',
'error': 'O(h⁴) global'
},
'Gaussian quadrature': {
'idea': 'Choose both nodes AND weights optimally',
'n-point rule': 'Exact for polynomials of degree ≤ 2n-1',
'Gauss-Legendre':'Nodes = roots of Legendre polynomials on [-1,1]',
'transform': '∫ₐᵇf(x)dx = (b-a)/2 ∫₋₁¹ f((b-a)t/2+(a+b)/2) dt',
'points_weights': {
'n=2': 'x=±1/√3, w=1',
'n=3': 'x=0,±√(3/5), w=8/9, 5/9'
}
},
'Adaptive quadrature': {
'idea': 'Subdivide intervals where error estimate is large',
'stopping': '|coarse - fine| < tolerance',
'Python': 'scipy.integrate.quad implements this'
},
'Monte Carlo integration': {
'formula': '∫f dx ≈ (b-a)/N Σf(xᵢ) where xᵢ uniform random',
'error': 'O(1/√N) regardless of dimension',
'best_for': 'High-dimensional integrals (d >> 1)'
}
}
def composite_simpson(f, a, b, n):
"""
Composite Simpson's rule: n must be even.
"""
if n % 2 != 0:
n += 1
h = (b - a) / n
x = np.linspace(a, b, n+1)
y = np.array([f(xi) for xi in x])
result = y[0] + y[-1]
result += 4 * np.sum(y[1:-1:2]) # odd indices
result += 2 * np.sum(y[2:-2:2]) # even indices (not endpoints)
return result * h / 3
def romberg_integration(f, a, b, max_rows=10):
"""
Romberg integration: Richardson extrapolation applied to trapezoidal rule.
R[n,m] = (4ᵐR[n,m-1] - R[n-1,m-1]) / (4ᵐ - 1)
"""
R = np.zeros((max_rows, max_rows))
# First column: trapezoidal with 2^n subintervals
for n in range(max_rows):
N = 2**n
h = (b-a)/N
x = np.linspace(a, b, N+1)
R[n,0] = h/2 * (f(a) + f(b) + 2*sum(f(xi) for xi in x[1:-1]))
# Richardson extrapolation
for m in range(1, max_rows):
for n in range(m, max_rows):
R[n,m] = (4**m * R[n,m-1] - R[n-1,m-1]) / (4**m - 1)
return R
def finite_differences():
return {
'Forward difference': '(f(x+h) - f(x))/h — O(h) error',
'Backward difference': '(f(x) - f(x-h))/h — O(h) error',
'Central difference': '(f(x+h) - f(x-h))/(2h) — O(h²) error',
'Second derivative': '(f(x+h) - 2f(x) + f(x-h))/h² — O(h²) error',
'Fourth order central': '(-f(x+2h)+8f(x+h)-8f(x-h)+f(x-2h))/(12h) — O(h⁴)',
'Optimal h': '~√ε_machine for first derivative ≈ 10⁻⁸',
'Trade-off': 'Smaller h: less truncation, more round-off'
}
def richardson_extrapolation(f, x, h, order=1):
"""
Richardson extrapolation to improve accuracy.
D(h) = f'(x) + C·hᵖ + O(h^(p+1))
D(h/2) = f'(x) + C·(h/2)ᵖ + ...
Extrapolated: f'(x) ≈ (2ᵖD(h/2) - D(h))/(2ᵖ-1)
"""
D_h = (f(x+h) - f(x-h)) / (2*h) # O(h²)
D_h2 = (f(x+h/2) - f(x-h/2)) / h # O(h²)
# Extrapolate to O(h⁴)
extrapolated = (4*D_h2 - D_h) / 3
return {
'D(h)': D_h,
'D(h/2)': D_h2,
'extrapolated': extrapolated
}
def numerical_linear_algebra():
return {
'LU decomposition': {
'factorization':'A = LU (L lower triangular, U upper triangular)',
'partial_pivot':'PA = LU for numerical stability',
'solve': 'Ax=b: Ly=Pb (forward), Ux=y (backward)',
'cost': 'O(n³) factorization, O(n²) per solve',
'use': 'Multiple right-hand sides with same A'
},
'Cholesky': {
'for': 'Symmetric positive definite A = LLᵀ',
'cost': 'Half the cost of LU',
'stability': 'Always stable for SPD matrices'
},
'QR decomposition': {
'factorization':'A = QR (Q orthogonal, R upper triangular)',
'methods': 'Gram-Schmidt, Householder (more stable), Givens',
'use': 'Least squares (more stable than normal equations)',
'eigenvalues': 'QR algorithm for all eigenvalues'
},
'Iterative methods (large sparse)': {
'Jacobi': 'xᵢ^(k+1) = (bᵢ - Σⱼ≠ᵢ aᵢⱼxⱼ^(k)) / aᵢᵢ',
'Gauss-Seidel': 'Use updated values immediately',
'SOR': 'Over-relaxation: ω∈(1,2) speeds Gauss-Seidel',
'CG': 'Conjugate Gradient for SPD (optimal Krylov)',
'GMRES': 'General non-symmetric systems',
'convergence': 'Depends on spectral radius ρ(iteration matrix)'
},
'Preconditioning': {
'idea': 'Solve M⁻¹Ax = M⁻¹b (better conditioned)',
'incomplete_LU':'ILU: sparse approximate factorization',
'diagonal': 'Jacobi preconditioner: M = diag(A)'
}
}
def power_iteration(A, num_iterations=100, tol=1e-10):
"""
Power iteration: finds largest eigenvalue and eigenvector.
v_{k+1} = Av_k / ||Av_k||
"""
n = A.shape[0]
v = np.random.rand(n)
v = v / np.linalg.norm(v)
eigenvalue = 0
for i in range(num_iterations):
Av = A @ v
eigenvalue_new = np.dot(v, Av)
v = Av / np.linalg.norm(Av)
if abs(eigenvalue_new - eigenvalue) < tol:
break
eigenvalue = eigenvalue_new
return eigenvalue, v
def ode_methods():
return {
'Euler (explicit)': {
'formula': 'y_{n+1} = y_n + h·f(t_n, y_n)',
'order': '1st order: error O(h)',
'stability': 'Conditionally stable: h < 2/|λ| for stiff',
'use': 'Simple demonstration, not production'
},
'Runge-Kutta 4': {
'k1': 'h·f(tₙ, yₙ)',
'k2': 'h·f(tₙ+h/2, yₙ+k1/2)',
'k3': 'h·f(tₙ+h/2, yₙ+k2/2)',
'k4': 'h·f(tₙ+h, yₙ+k3)',
'formula': 'y_{n+1} = y_n + (k1+2k2+2k3+k4)/6',
'order': '4th order: error O(h⁴)',
'use': 'Standard for non-stiff problems'
},
'Dormand-Prince (RK45)': {
'idea': 'Embedded 4th and 5th order for error estimation',
'adaptive': 'Adjust h based on ||y₅-y₄||',
'implementation':'scipy.integrate.solve_ivp default'
},
'Stiff equations': {
'problem': 'Solution varies on multiple timescales',
'symptom': 'Explicit methods require tiny h',
'solution': 'Implicit methods (backward Euler, trapezoidal, BDF)',
'backward_euler':'y_{n+1} = y_n + h·f(t_{n+1}, y_{n+1}) (solve nonlinear)',
'BDF': 'Backward Differentiation Formulas (scipy ode15s equivalent)'
},
'Multistep methods': {
'Adams-Bashforth':'Explicit, uses previous steps',
'Adams-Moulton': 'Implicit, higher accuracy',
'BDF': 'Stiff solvers, variable order'
}
}
def rk4_solver(f, y0, t_span, h):
"""
Classic RK4 ODE solver.
"""
t0, tf = t_span
t = np.arange(t0, tf+h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
for i in range(len(t)-1):
k1 = h * np.array(f(t[i], y[i]))
k2 = h * np.array(f(t[i]+h/2, y[i]+k1/2))
k3 = h * np.array(f(t[i]+h/2, y[i]+k2/2))
k4 = h * np.array(f(t[i]+h, y[i]+k3))
y[i+1] = y[i] + (k1 + 2*k2 + 2*k3 + k4) / 6
return t, y
def pde_methods():
return {
'Finite Difference Method (FDM)': {
'idea': 'Replace derivatives with finite differences on grid',
'heat_explicit':'u^(n+1)_i = u^n_i + r(u^n_{i+1}-2u^n_i+u^n_{i-1}), r=αΔt/Δx²',
'stability': 'Explicit heat: r ≤ 1/2 (CFL condition)',
'implicit': 'Always stable (Crank-Nicolson, backward Euler)',
'wave': 'CFL: c·Δt/Δx ≤ 1 for stability'
},
'Finite Element Method (FEM)': {
'idea': 'Variational formulation, piecewise polynomial basis',
'weak_form': 'Multiply PDE by test function, integrate by parts',
'stiffness': 'Assemble global stiffness matrix K and load vector f',
'solve': 'Ku = f (linear system)',
'advantage': 'Complex geometry, natural boundary conditions',
'software': 'FEniCS, COMSOL, ANSYS, Abaqus'
},
'Finite Volume Method (FVM)': {
'idea': 'Integrate PDE over control volumes (conservation form)',
'advantage': 'Exactly conservative (important for fluids)',
'software': 'OpenFOAM, SU2'
},
'Spectral methods': {
'idea': 'Represent solution as global basis functions (Fourier, Chebyshev)',
'advantage': 'Exponential convergence for smooth solutions',
'use': 'Smooth periodic or simple geometry problems'
}
}
def optimization_methods():
return {
'Unconstrained optimization': {
'Gradient descent': 'x_{k+1} = x_k - α∇f(x_k)',
'Newton\'s method': 'x_{k+1} = x_k - [∇²f(x_k)]⁻¹∇f(x_k)',
'BFGS': 'Quasi-Newton: approximate Hessian update',
'L-BFGS': 'Limited memory BFGS (large scale)',
'Conjugate gradient': 'For quadratic problems and extensions',
'Line search': 'Choose α to satisfy Wolfe conditions'
},
'Constrained optimization': {
'Lagrange multipliers': '∇f = λ∇g (equality constraints)',
'KKT conditions': 'First-order necessary for inequality constraints',
'Penalty methods': 'Add constraint violation to objective',
'Interior point': 'Stay strictly feasible, log barrier',
'Active set': 'Identify which constraints are active'
},
'Global optimization': {
'Simulated annealing': 'Accept worse solutions with decreasing probability',
'Genetic algorithms': 'Evolutionary: selection, crossover, mutation',
'Particle swarm': 'Swarm intelligence',
'Basin hopping': 'Random perturbation + local optimization'
}
}
def numerical_software():
return {
'Python': {
'NumPy': 'Array operations, linear algebra',
'SciPy': 'Root finding, integration, ODEs, optimization, linear algebra',
'Matplotlib': 'Visualization',
'SymPy': 'Symbolic mathematics'
},
'Key SciPy functions': {
'scipy.optimize.brentq': 'Root finding (robust)',
'scipy.optimize.fsolve': 'System of nonlinear equations',
'scipy.integrate.quad': 'Adaptive quadrature',
'scipy.integrate.solve_ivp':'ODE solver (RK45 default)',
'scipy.linalg.solve': 'Linear system (uses LAPACK)',
'scipy.linalg.lstsq': 'Least squares',
'scipy.linalg.eig': 'Eigenvalues',
'scipy.sparse.linalg.spsolve': 'Sparse linear system'
},
'MATLAB equivalents': {
'ode45': 'Non-stiff ODE (Dormand-Prince RK45)',
'ode15s': 'Stiff ODE (BDF/NDF)',
'fzero': 'Root finding',
'fminunc': 'Unconstrained optimization',
'linprog': 'Linear programming',
'\\': 'Linear system solve (backslash operator)'
}
}
| Pitfall | Fix |
|---|---|
| Small h always better for differentiation | Round-off increases as h→0; use h≈√ε for first derivative |
| Newton always converges | Needs good initial guess near root; can diverge or cycle |
| Large n polynomial interpolation | Use splines or Chebyshev nodes to avoid Runge's phenomenon |
| Explicit method for stiff ODE | Use implicit method (BDF, Crank-Nicolson) |
| Normal equations for least squares | Use QR decomposition; normal equations square the condition number |
| Ignoring condition number | κ ≈ 10⁸ means losing ~8 digits of accuracy |
まだレビューはありません。使ってみた感想をお寄せください。
概要と使いどころ
Expert-level abstract algebra knowledge. Use when working with groups, rings, fields, modules, Galois theory, representation theory, or algebraic structures. Also use when the user mentions 'group', 'ring', 'field', 'homomorphism', 'isomorphism', 'normal subgroup', 'quotient group', 'Sylow theorem', 'Galois theory', 'polynomial ring', 'ideal', 'module', 'vector space', or 'representation theory'.
日本語の概要は準備中です。原文の説明を表示しています。
Expert-level acoustics knowledge. Use when working with sound waves, acoustic propagation, resonance, standing waves, Doppler effect, room acoustics, ultrasound, noise control, musical acoustics, or architectural acoustics. Also use when the user mentions 'sound wave', 'decibel', 'frequency', 'resonance', 'standing wave', 'Doppler effect', 'acoustic impedance', 'reverberation', 'ultrasound', 'noise control', 'harmonics', or 'sound intensity'.
日本語の概要は準備中です。原文の説明を表示しています。
Expert-level aerodynamics covering subsonic and supersonic flow, lift and drag, airfoil theory, boundary layers, compressible flow, and CFD methods.
日本語の概要は準備中です。原文の説明を表示しています。
Expert-level air quality engineering covering pollutant sources, atmospheric dispersion, air quality standards, emission controls, and air quality modeling.
日本語の概要は準備中です。原文の説明を表示しています。
Expert-level aircraft structures covering structural design philosophy, wing and fuselage analysis, fatigue and damage tolerance, composite structures, and certification.
日本語の概要は準備中です。原文の説明を表示しています。
Expert-level algorithms and complexity theory covering sorting, searching, graph algorithms, dynamic programming, greedy algorithms, NP-completeness, and algorithm design paradigms.
日本語の概要は準備中です。原文の説明を表示しています。