Найти решение трансцендентного уравнения с помощью решающего кода - PullRequest
0 голосов
/ 05 декабря 2018

я пытаюсь найти решение трансцендентного уравнения вида:

(arccos (tan (lat) / tan (inc))) - ((n * pi) + (arccos (sin (lat) / sin (inc)))) * 1 / q = 0

(Это часть статьи по вычислению пересмотра орбиты ...) Я знаю, что параметры lat, n и qможет быть предоставлено, и для параметра inc должно быть найдено только решение.

Не всегда есть решение ... n и q должны быть разрешены некоторые значения, но, например, если lat равно -60,тогда n = 1 и 3/8 <1 / q <.5 должны работать, чтобы иметь решение inc между abs (lat) и 180-abs (lat), например, если lat имеет значение -60 градусов, то incдолжно иметь значение между 60 и 120 градусами, это то, что говорят физики ... </p>

Я попробовал следующий код, основанный на предложенных решениях для других трансцендентных уравнений здесь в stackoverflow, но я не могу получить значение inc в пределахожидаемые значения, и в большинстве случаев ошибка заключается в том, что arcos (tan (lat) / tan (inc)) не является ожидаемым значением междуn -1 и 1, но когда это не так, решатель все еще не может решить и предоставить ожидаемое значение для inc.

Я использую этот код для тестирования:

from scipy.optimize import minimize
from numpy import arccos, sin, tan, radians,pi,degrees
def opt_fun(n,lat,q,inc):
    if abs(tan(lat)/tan(inc))<1.000:
        print(lat,inc)
        f=(arccos(tan(lat)/tan(inc)))-((n*pi)+(arccos(sin(lat)/sin(inc))))*1/q
    else:
        f=0
    return(f)

n=1
lat=radians(-60.0)
q = 2
for i in range(1,100):
    q+=.1
    try:
        res = minimize(lambda inc: opt_fun(n,lat,q,inc), x0=abs(lat))

        # Check if the optimization was successful
        print(res.success)
        if res.success:
            # expected >> True

            # Extract the root from the minimization result
            print(degrees(res.x[:]))
            # expected: abs(lat) < inc < 180-abs(lat)
    except:
        pass

Если вызапустив его, вы увидите, что единственное возможное решение - это когда abs (tan (lat) / tan (inc)) <1.000, но значения inc не находятся в пределах 0,2PI, как я ожидаю для угла в радианах.</p>

Должен ли я ограничить отклонение x0 между радианами (абс (лат)) и радианами (180 лат)?как?или есть какая-то другая проблема, которую нужно решить?

Если мы подавляем условие, которое не позволяет входным данным arccos быть выше -1,1, то появляется следующая ошибка: «RuntimeWarning: недопустимое значение, встречающееся в arccos f= (arccos (tan (lat) / tan (inc))) - ((n * pi) + (arccos (sin (lat) / sin (inc)))) * 1 / q "

IЯ не очень беспокоюсь об эффективности этого проекта кода, но мне нужны хотя бы некоторые результаты, чтобы продолжить тестирование.Заранее спасибо.С уважением!

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...