Матлаб: реализовать функцию для узлов произвольной решетки

[GAA]

А класс "рациональные" есть в Matlab? Если же использовать символьные, то — нет, см. справку/документацию (можно в сети sym/mod). Или я что-то не понял?

[maximkarimov]

arseniiv в сообщении #1514827 писал(а):

А когда решётка квадратная, шаг не важен

Да, решетка всегда квадратная и всегда обычная (без вращения вокруг начала координат).
Рассмотрим конкретный, простейший случай:
1) луч $X$ выходит из начала координат и проходит через точку А с координатами ($\cos(2),\sin(2)$);
2) решетка имеет постоянный шаг $\sqrt{2}$.
Требуется найти координаты узла, через который проходит $X$ и расстояние от которого до начала координат минимально.

-- 17.04.2021, 22:36 --

GAA в сообщении #1514828 писал(а):

Второй вариант заведомо нерабочий

Да вроде работает:

Используется синтаксис Matlab M

>> x=sym(10*pi)
 
x =
 
10*pi
 
>> y=sym(5*pi)
 
y =
 
5*pi
 
>> isknot(x,y,pi)

ans =

     1

[Dmitriy40]

А зачем здесь иррациональности? Берёте угол наклона луча в единицах решётки $\frac{\sin(2)/\sqrt{2}}{\cos(2)/\sqrt{2}}$ и смотрите через какую целочисленную точку $\frac{y}{x}$ он теперь проходит. Целые числа нумеруют шаги решётки.
Если шаг решётки разный по разным осям, то $\sqrt{2}$ в числителе и знаменателе будут разными, по соответствующей оси.
И если угол наклона получился иррациональным, то подходящих узлов решётки нет, с каким бы постоянным шагом по каждой оси она не была. Просто потому что целочисленная решётка может содержать лишь рациональные углы.
Целый, рациональный или иррациональный шаг решётки роли не играет, важно чтобы угол наклона луча был рациональным в шагах решётки.
В принципе даже не обязательно именно постоянный шаг решётки по осям, можно и переменный, но тогда функции задающие решётку $x=f_x(n), y=f_y(n)$ в простейшем случае должны быть одинаковы для всех осей. Или точнее такими, чтобы угол наклона луча получался рациональным. Вот Вам и условие на параметры решётки.

PS. Конечно везде выше вместо "угол" надо бы писать "тангенс угла" ...

[maximkarimov]

Dmitriy40 в сообщении #1514837 писал(а):

Рациональный, целый или иррациональный шаг решётки роли не играет, важно чтобы угол наклона луча был рациональным в шагах решётки.

Мм.. на вскидку не соображу как идентифицировать. Ну хорошо, это необходимое условие на прохождение луча хотя бы через один узел решетки, так? Но ведь нам нужно найти координаты конкретного узла...

[Dmitriy40]

maximkarimov
Как только получите конкретное рациональное число, то превратить его в координаты тривиально: числитель это номер линии решётки по оси $y$, знаменатель это номер линии решётки по оси $x$. Пересчитать в координаты легко, просто домножить на шаг решётки по каждой оси. И да, это ближайший к началу координат узел, просто по свойствам рациональных чисел (дробь несократима).
А насчёт идентифицировать ... Видимо надо проследить чтобы все иррациональности сократились. Хотя бы при делении каждой из координат на величину шага решётки по каждой оси.

-- 17.04.2021, 22:09 --

Например, луч через точку $(\ln(10000);10\sin(2))$ пройдёт и через узел решётки $(2;5)$ (это номера линий решётки, не координаты) при шагах решётки по осям $(\ln(10);\sin(2))$. Потому что $\frac{\ln(10000)/\ln(10)}{10\sin(2)/\sin(2)}=\frac{4}{10}=\frac{2}{5}$, т.е. рационально. Тут шаг постоянный, но не одинаковый по осям.

[GAA]

maximkarimov в сообщении #1514832 писал(а):

Да вроде работает:

В документации пример такого случая не приведен. В старых версиях, например Matlab 2013, такой вариант не работал. Но в 4 версии Оtave такой вариант работает. Видимо допилили [в Matlab] и у Вас достаточно свежая версия. К рассматриваемой Вами задаче этот вопрос, конечно не относится, но момент для меня интересный. Буду знать, спасибо.

[maximkarimov]

GAA в сообщении #1514844 писал(а):

К рассматриваемой Вами задаче этот вопрос, конечно не относится, но момент для меня интересный. Буду знать, спасибо.

Вам спасибо!

-- 17.04.2021, 23:25 --

Dmitriy40 в сообщении #1514841 писал(а):

maximkarimov
Например, луч через точку $(\ln(10000);10\sin(2))$ пройдёт и через узел решётки $(2;5)$ (это номера линий решётки, не координаты) при шагах решётки по осям $(\ln(10);\sin(2))$. Потому что $\frac{\ln(10000)/\ln(10)}{10\sin(2)/\sin(2)}=\frac{4}{10}=\frac{2}{5}$, т.е. рационально. Тут шаг постоянный, но не одинаковый по осям.

Кажись начинаю понимать... Спасибо!
Забегаю вперед, но все же спрошу - а как быть в случае, если решетка задается частичными суммами (сходящегося) ряда?

[GAA]

По поводу задачи с лучом. Как и писал выше arseniiv, боюсь в вопросе точного решения Matlab не помощник. Можно поставить вопрос о прохождении вблизи некоторого узла. В Matlab есть функция rat, которая подбирает приближенное значение дроби. В моей версии rat для символьных нет, для чисел с плавающей точкой можно попробовать:

Используется синтаксис Matlab M

[n, d] = rat(cos(2)/sin(2), 1e-15)
n = -37774882
d = 82539623

Видно, что, по крайней мере, вблизи начала координат луч через узел не проходит.
В свежих версиях Matlab rat для символьных должна быть. Там можно получить огромные числа в числителе и знаменателе.

Если шаг не постоянный (хоть ряды, хоть что-то другое), то, как и писал выше arseniiv, прямой перебор до достаточно больших номеров узлов на базе символьной арифметики.

[GAA]

GAA в сообщении #1514844 писал(а):

Но в 4 версии Оtave такой вариант работает.

Естественно без double

Используется синтаксис Matlab M

function [res] = isknot(x,s)
if mod(x,s)==0
    res=1;
else
    res=0;
end

Используется синтаксис Matlab M

>> x = sym('(1+1/10000000000000000000)*pi')
x = (sym)

  10000000000000000001*pi
  -----------------------
    10000000000000000000

>> isknot(x, sym(pi))
ans = 0
>> x = sym('10*pi');
>> isknot(x, sym(pi))
ans =  1

Необходимость использовать double в свежих версиях Matlab выглядит странно: нет синтаксической необходимости и уменьшает точность вычислений.

[Dmitriy40]

maximkarimov в сообщении #1514845 писал(а):

Забегаю вперед, но все же спрошу - а как быть в случае, если решетка задается частичными суммами (сходящегося) ряда?

Я не очень понимаю при чём здесь ряд, но если у Вас каждый шаг решётки является некоей функцией от его номера, типа $f_{x,y}(n)$, то скорее всего придётся потребовать чтобы отношение шагов по осям $f_y(n)/f_x(n)$ было одинаковым для всех $n$, сами шаги могут быть любыми ненулевыми. Потому что это будет просто масштабирование с одинаковым коэффициентом по осям на каждом шаге, а оно углы не меняет. Но тогда отпадает возможность использовать иррациональные отношения шагов по осям.
Возможно в тензорной форме и можно что-то придумать (есть же криволинейные координаты или распространение луча в оптически неоднородной среде), но тут я полный нуль.

[arseniiv]

GAA в сообщении #1514831 писал(а):

А класс "рациональные" есть в Matlab?

Не знаю, не специалист по этой СКА, но по идее это должно бы быть во всех достаточно общепрофильных (в Mathematica есть, плюс например язык Python умеет это сам по себе без дополнительных пакетов), так как карман не тянет.

[GAA]

Octave/Matlab как бы не СКА, а Maple — СКА. В Octave 4 нет в поставке класса рациональные, но для символьных

Используется синтаксис Matlab M

>> mod(sym('2+1/2'), sym('1/3'))
ans = (sym) 1/6

Тогда как в Maple 15 обещают в справке только для целых модулей:

> (2+1/2) mod 1/2;
Error, invalid input: modp received 1/2, which is not valid for its 2nd argument, m
> `mod`(2+1/2,  1/2);
Error, invalid input: modp received 1/2, which is not valid for its 2nd argument, m

Если значения входных аргументов иррациональные, то возвращаемые результаты в Octave 4, за исключением очень простых случаев, несколько неожиданны

 >> mod(sym('(2+1/2)*pi'), sym('pi'))
ans = (sym)
  pi
  --
  2
>> mod(sym('2*pi'), sym('pi'))
ans = (sym) 0
>> mod(sym('3*cos(3)'), sym('sin(2)'))
ans = (sym) 3*cos(3) + 4*sin(2)

[Octave, грубо говоря, использует SymPy для символьных вычислений.] На удивление в некоторых вырожденных случаях в Maple работает для иррациональных модулей.

> 2*Pi mod Pi;
                   0
> # but!
> (2+1/2)*Pi mod Pi;
               modp(5/2*Pi,Pi)

[В Мaple удобно, если не получилось вычислить, получить исходное выражение.]
В общем, нужны ли хоть какие-то случаи с иррациональными вопрос технический, при программировании что-то может быть полезным.

[maximkarimov]

GAA в сообщении #1514851 писал(а):

боюсь в вопросе точного решения Matlab не помощник. ... Видно, что, по крайней мере, вблизи начала координат луч через узел не проходит.
В свежих версиях Matlab rat для символьных должна быть. Там можно получить огромные числа в числителе и знаменателе

Вот что выдала моя версия Матлаб (2017b):

Используется синтаксис Matlab M

>> rat(cos(2)/sin(2))

ans =

    '-0 + 1/(-2 + 1/(-5 + 1/(-2 + 1/(-2 + 1/(-9)))))'

>> sym(-0 + 1/(-2 + 1/(-5 + 1/(-2 + 1/(-2 + 1/(-9))))))
 
ans =
 
-254/555

>> a=sym(2)
 
a =
 
2

>> x=-254*sqrt(a)

x =
 
-254*2^(1/2)
 
>> y=255*sqrt(a)
 
y =
 
255*2^(1/2)
 
>> isknot(x,y,sqrt(a))

ans =

     1

 

Ответ: луч X проходит через узел решетки с координатами (-254$\sqrt{2}$, 255$\sqrt{2}$), который является ближайшим к началу координат. Правильно?

[GAA]

Вы же в справке можете посмотреть погрешность (по умолчанию) с которой rat приближает иррациональное выражение. Поэтому "не проходит через узел такой-то", а "проходит вблизи узла такого-то на расстоянии не большем чем …" :)

[maximkarimov]

Dmitriy40 в сообщении #1514837 писал(а):

И если угол наклона получился иррациональным, то подходящих узлов решётки нет, с каким бы постоянным шагом по каждой оси она не была. Просто потому что целочисленная решётка может содержать лишь рациональные углы.
Целый, рациональный или иррациональный шаг решётки роли не играет, важно чтобы угол наклона луча был рациональным в шагах решётки.
В принципе даже не обязательно именно постоянный шаг решётки по осям, можно и переменный, но тогда функции задающие решётку $x=f_x(n), y=f_y(n)$ в простейшем случае должны быть одинаковы для всех осей. Или точнее такими, чтобы угол наклона луча получался рациональным. Вот Вам и условие на параметры решётки.

PS. Конечно везде выше вместо "угол" надо бы писать "тангенс угла" ...

Сорри, что до меня как до жирафа доходит!) Спасибо большое! С этим типом решеток "покончено", приступаю к формулировке условий на решетки другого типа. Напишу в этой теме - вдруг там тоже без перебора узлов что-то можно будет понять...