Я пытаюсь минимизировать действие, пытаясь выяснить, как ведут себя параметры. Мое действие - это , и я стараюсь, чтобы оно стало этим для численной оценки. Различая это, я нахожу уравнения для dS / dA_m = 0 и dS / dphi_m = 0 , и я хочу использовать метод Ньютона-Рафсона, чтобы обнаружить функции A и phi, зависящие от r которые решают эти уравнения.
Я определил функции для моей установки, в том числе нахождение якобиана для метода Ньютона, следующим образом:
N = 5
dr = 0.01
A = np.arange(1,N+1,1)
phi = np.arange(1,N+1,1)
def V(phi,phi_M):
return -3-3/2*phi**2-1/3*phi**4+(1/(3*phi_M**2)+1/(2*phi_M**4))*phi**6-1/(12*phi_M**4)*phi**8
def dV(phi,phi_M):
return -3*phi-4/3*phi**3+(2/(phi_M**2)+3/(phi_M**4))*phi**5-2/(3*phi_M**4)*phi**7
def ddV(phi,phi_M):
return -3-4*phi**2+5*(2/(phi_M**2)+3/(phi_M**4))*phi**4-14/(3*phi_M**4)*phi**6
def dS_dphi(A,phi,phi_M,dr):
vector = []
for i in range(0,N-1):
vector.append(-dr*((phi[i]-phi[i-1])*np.exp(4*A[i-1])/dr**2-(phi[i+1]-phi[i])*np.exp(4*A[i])/dr**2+dV(phi[i],phi_M)*np.exp(4*A[i])))
return vector
def dS_dA(A,phi,phi_M,dr):
vector = []
for i in range(0,N-1):
vector.append(dr*(4*A[i]*np.exp(4*A[i])*(3*(A[i+1]-A[i])**2/dr**2-(phi[i+1]-phi[i])**2/(2*dr**2)-V(phi[i],phi_M))
+6*(np.exp(4*A[i-1])*(A[i]-A[i-1])/dr**2-np.exp(4*A[i])*(A[i+1]-A[i])/dr**2)))
return vector
def d2S_dphildphim(A,phi,phi_M,dr):
matrix = []
for i in range(0,N-1):
row = []
for j in range (0,N-1):
if i == j:
row.append(-dr*(np.exp(4*A[i-1])/dr**2+np.exp(4*A[i])/dr**2+ddV(phi[i],phi_M)*np.exp(4*A[i])))
elif i == j+1:
row.append(dr*(np.exp(4*A[i-1])/dr**2))
elif i == j-1:
row.append(dr*(np.exp(4*A[i])/dr**2))
else:
row.append(0)
matrix.append(row)
return matrix
def d2S_dAldphim(A,phi,phi_M,dr):
matrix = []
for i in range(0,N-1):
row = []
for j in range (0,N-1):
if i == j:
row.append(dr*(4*(phi[i+1]-phi[i])/dr**2*np.exp(4*A[i])-4*np.exp(4*A[i])*dV(phi[i],phi_M)))
elif i == j-1:
row.append(-dr*(4*(phi[i+1]-phi[i])*np.exp(4*A[i])/dr**2))
else:
row.append(0)
matrix.append(row)
return matrix
def d2S_dphildAm(A,phi,phi_M,dr):
matrix = []
for i in range(0,N-1):
row = []
for j in range (0,N-1):
if i == j:
row.append(dr*4*np.exp(4*A[i])*((phi[i+1]-phi[i])/dr**2-dV(phi[i],phi_M)))
elif i == j+1:
row.append(-dr*4*np.exp(4*A[i-1])*(phi[i]-phi[i-1])/dr**2)
else:
row.append(0)
matrix.append(row)
return matrix
def d2S_dAldAm(A,phi,phi_M,dr):
matrix = []
for i in range(0,N-1):
row = []
for j in range (0,N-1):
if i == j:
row.append(dr*(16*np.exp(4*A[i])*(3*(A[i+1]-A[i])**2/dr**2-(phi[i+1]-phi[i])**2/(2*dr**2)-V(phi[i],phi_M)
-3*(A[i+1]-A[i])/dr**2+3/(8*dr**2))+6*np.exp(4*A[i-1])/dr**2))
elif i == j+1:
row.append(dr*6*np.exp(4*A[i-1])*(4*(A[i]-A[i-1])/dr**2-1/dr**2))
elif i == j-1:
row.append(dr*6*np.exp(4*A[i])*(4*(A[i+1]-A[i])/dr**2-1/dr**2))
else:
row.append(0)
matrix.append(row)
return matrix
def Jacobian(A,phi,phi_M,dr):
top = np.concatenate((np.array(d2S_dAldAm(A,phi,phi_M,dr)),np.array(d2S_dphildAm(A,phi,phi_M,dr))),axis=1)
bottom = np.concatenate((np.array(d2S_dAldphim(A,phi,phi_M,dr)),np.array(d2S_dphildphim(A,phi,phi_M,dr))),axis=1)
return np.concatenate((top,bottom),axis=0)
, но я не совсем уверен, как реализовать их для фактического их решения. для функций. Это то, что я взял откуда-то еще и попытался преобразовать в свою форму, которая, как я знаю, неверна:
x = -2
y = -1.5
i1 = 0
while i1<10:
F= np.concatenate((np.array(dS_dphi(A,phi,1,dr)),np.array(dS_dA(A,phi,1,dr))),axis=0)
theta = np.sum(F)
J = Jacobian(A,phi,1,dr)
Jinv = inv(J)
xn = np.array([[x],[y]])
xn_1 = xn - (Jinv*F)
x = xn_1.item(0)
y = xn_1.item(1)
#~ print theta
print(xn)
i1 = i1+1
Может ли кто-нибудь указать мне направление или решение, которое поможет мне решить эту проблему? ? Большое спасибо! :)