{}
from math import sqrt
def J_prime(x, m, n, ax, ay, bx, by):
""" Derivada de J(x) con respecto a x """
DA = sqrt((x - ax) ** 2 + (m * x + n - ay) ** 2)
DB = sqrt((x - bx) ** 2 + (m * x + n - by) ** 2)
return ((x - ax) + m * (m * x + n - ay)) / DA + ((x - bx) + m * (m * x + n - by)) / DB
def J_double_prime(x, m, n, ax, ay, bx, by):
""" Segunda derivada de J(x) con respecto a x """
y = m * x + n
dy_dx = m
# Compute terms for PA
dx_pa = x - ax
dy_pa = y - ay
DA = sqrt(dx_pa**2 + dy_pa**2)
numPA = dx_pa + m * dy_pa # (x - ax) + m*(y - ay)
term1_numerator = (1 + m**2) * DA**2 - numPA**2
term1 = term1_numerator / DA**3
# Compute terms for PB
dx_pb = x - bx
dy_pb = y - by
DB = sqrt(dx_pb**2 + dy_pb**2)
numPB = dx_pb + m * dy_pb # (x - bx) + m*(y - by)
term2_numerator = (1 + m**2) * DB**2 - numPB**2
term2 = term2_numerator / DB**3
return term1 + term2
def find_root_bisection(m, n, ax, ay, bx, by, a, b, tol=1e-6, max_iter=50):
""" Método de bisección para encontrar un intervalo donde J'(x) cambia de signo """
fa = J_prime(a, m, n, ax, ay, bx, by)
fb = J_prime(b, m, n, ax, ay, bx, by)
if fa * fb > 0:
print("Bisección falló: No hay cambio de signo en el intervalo.")
return None
for _ in range(max_iter):
c = (a + b) / 2
fc = J_prime(c, m, n, ax, ay, bx, by)
if abs(fc) < tol:
return c
if fa * fc < 0:
b, fb = c, fc
else:
a, fa = c, fc
return (a + b) / 2
def newton_raphson(x0, m, n, ax, ay, bx, by, tol=1e-6, max_iter=10):
""" Método de Newton-Raphson con control de estabilidad """
x = x0
for i in range(max_iter):
fx = J_prime(x, m, n, ax, ay, bx, by)
dfx = J_double_prime(x, m, n, ax, ay, bx, by)
if abs(dfx) < tol:
print("Segunda derivada cercana a cero. Newton-Raphson no es estable.")
return x # Devuelve el mejor valor actual sin continuar
x_new = x - fx / dfx
if abs(x_new - x) < tol:
break
x = x_new
return x
# Parámetros de la recta y los puntos
m, n = 0.5, 1
ax, ay = 1, 3
bx, by = 4, 5
intervalo = (1, 5) # Rango donde buscamos la raíz
# Paso 1: Buscar una raíz con bisección
x0 = find_root_bisection(m, n, ax, ay, bx, by, *intervalo)
if x0 is not None:
print(f"Raíz inicial encontrada por bisección: x0 = {x0:.6f}")
# Paso 2: Refinar con Newton-Raphson
x_final = newton_raphson(x0, m, n, ax, ay, bx, by)
y_final = m * x_final + n
print(f"Solución final: x = {x_final:.6f}, y = {y_final:.6f}")
else:
print("No se encontró raíz en el intervalo dado.")