Магнитное поле однородно намагниченного цилиндра, как проверить формулу из ПБМ2?

peregoudov · Сообщение **peregoudov** » 28 ноя 2018, 14:53

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

Используя верную вне магнита формулу

${\bf B}=({\bf m}\nabla)\nabla\psi,\quad \psi({\bf r})=\int\frac{d^3r'}{|{\bf r}-{\bf r}'|}$

и разложение фундаментального решения по функциям Бесселя

$\frac1{|{\bf r}-{\bf r}'|}= \sum_{m=-\infty}^{+\infty}e^{im(\phi-\phi')} \int_0^\infty J_m(kr)J_m(kr')e^{-k|z-z'|}\,dk$

(в правой части $(r,\phi,z)$ --- цилиндрические координаты), дошел вот до какого выражения

${\bf B}(r,\phi,z)= \int_0^\infty R\,dk\,g(kz,kh)J_1(kR)\biggl[ (-{\bf e}_rm_x\cos\phi+{\bf e}_zm_z)J_0(kr)+{}$
${}+({\bf e}_zm_x\cos\phi+{\bf e}_rm_z)J_1(kr)+ ({\bf e}_r\cos\phi+{\bf e}_\phi\sin\phi)m_x \frac{J_1(kr)}{kr}\biggr].$

Здесь радиус цилиндра обозначен R, высота h, $m_y=0$

,

$g(kz,kh)=\int_{-h/2}^{h/2}e^{-k|z-z'|}\,k\,dz'$

--- сводится к набору экспонент.

Видно, что возникает три эталонных интеграла, первый из них

$\int_0^\infty dx\,e^{-ax}J_1(bx)J_0(cx).$

В Прудникове---Брычкове---Маричеве (том 2, раздел 2.12.38) приводится такой ответ

$-\frac{ak}{2\pi b^{3/2}c^{1/2}}\left[ {\bf K}(k)\pm\frac1{2b}\Lambda_0(\phi,\psi)\right] +\left\{\begin{matrix}0\\1/b\end{matrix}\right\},$

где

$k=\sin\phi=\frac{2\sqrt{bc}}{\sqrt{a^2+(b+c)^2}},\quad \sin\psi=\frac a{\sqrt{a^2+(b-c)^2}},$

верхний знак (и верхний член в фигурных скобках) для $b<c$

, нижний --- для $b>c$

. Определение функции $\Lambda_0$ я еле нашел на стр. 662

$\Lambda_0(\phi,k)=\frac2\pi\left[ {\bf E}(k)F(\phi,\sqrt{1-k^2})+ {\bf K}(k)E(\phi,\sqrt{1-k^2})- {\bf K}(k)F(\phi,\sqrt{1-k^2})\right],$

причем на стр. 657 и 659 расшифрованы стандартные эллиптические интегралы

$\{E,F\}(\phi,k)= \int_0^\phi(1-k^2\sin^2t)^{\pm1/2}\,dt,\quad \{{\bf E},{\bf K}\}(k)=\{E,F\}(\pi/2,k).$

А вот дальше проблема: при прямом подсчете интеграла для (a,b,c)=(1,1,2) получаем 0.05185272, а при расчете по формуле из ПБМ получаем -0.266674632. При этом другая формула, для интеграла с двумя $J_1$

, правильная. В чем тут может быть дело?

P. S. Тему можно в некотором смысле считать продолжением этой.

peregoudov · Сообщение **peregoudov** » 10 дек 2018, 13:52

Я, тем временем, обнаружил, что Maple, похоже, неверно преобразует эллиптические интегралы. Потому что на вопрос посчитать $I(1,1,2)$

, где

$I(a,b,c)=\int_{-1}^1\frac{a-ict}{\sqrt{1-t^2}\,\sqrt{(a-ict)^2+b^2}}\,dt,$

он отвечает значением $-0.911685236-1.997134816\,i$ , тогда как правильный ответ $2.978692530$

.

Взялся я считать сам. Прежде всего представляем $J_0$

в виде интерала

$J_0(x)=\frac1{2\pi}\oint\frac{e^{ixt}\,dt}{\sqrt{1-t^2}}.$

Тогда

$I_{1,0}=\int_0^\infty e^{-ax}J_1(bx)J_0(cx)\,dx= \frac1{2\pi}\oint\frac{dt}{\sqrt{1-t^2}}\int_0^\infty e^{-(a-ict)x}J_1(bx)\,dx.$

Теперь пользуемся лапласовским образом

$\int_0^\infty e^{-px}J_1(x)\,dx=1-\frac p{\sqrt{p^2+1}},$

в результате чего и приходим к эллиптическому интегралу

$I_{1,0}=\frac1b-\frac1{2\pi b}\oint\frac{a-ict}{\sqrt{1-t^2}\,\sqrt{(a-ict)^2+b^2}}\,dt.$

Теперь задача в том, чтобы привести это к стандартным эллиптическим интегралам. Подынтегральное выражение имеет точки ветвления $t=\pm1$ , $\pm b/c-ia/c$ , лежащие, в силу симметрии, на одной окружности. Обозначая нижнюю точку пересечения этой окружности с мнимой осью через $-id$

, сконструируем дробно-линейное преобразование $t\to z$ , которое оставляло бы $\pm1$ неподвижными, а $-id$

переводило бы в бесконечность. Тогда указанная окружность на плоскости переменной $z$

будет вещественной прямой, а интеграл приведется к стандартному виду. Указанное преобразование выписывается элементарно

$t=\frac{1-idz}{z-id}.$

Величину $d$

проще всего найти, подставив преобразование в интеграл и потребовав, чтобы под корнем были только четные степени $z$

$d=\frac{a^2+b^2-c^2+\sqrt{(a^2+b^2-c^2)^2+4a^2c^2}}{2ac}.$

Сам же интеграл принимает вид

$I_{1,0}=\frac1b-\frac1{2\pi b}\oint\frac{\sqrt{1+d^2}\,[(-a+cd)z+i(ad+c)]\,dz} {(-z+id)\sqrt{1-z^2}\,\sqrt{(a^2+b^2)d^2+c^2+2acd-(a^2+b^2+c^2d^2-2acd)z^2}}$

и выражается через стандартные интегралы

$I_{1,0}=\frac1b-\frac{2\sqrt{1+d^2}}{\pi b\sqrt{(a^2+b^2)d^2+c^2+2acd}}\left[ (a-cd)K(k)+\frac{c(1+d^2)}d\Pi(-1/d^2,k)\right],$

здесь

$k=\sqrt{\frac{a^2+b^2+c^2d^2-2acd}{(a^2+b^2)d^2+c^2+2acd}}.$

И вот это выражение уже дает правильное значение на Maple.

Только я теперь ума не приложу, как Прудникову и Ко удалось выразить интеграл третьего рода через первый и второй... Или есть возможность как-то по другому провести преобразование интеграла?

peregoudov · Сообщение **peregoudov** » 11 дек 2018, 15:17

Ничего не понимаю... Читаю русскую Википедию, и вижу формулы, очень похожие на Прудникова и Ко. Но в английской версии страницы этих формул нет!

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

P. S. Обнаружил, откуда все списано: Абрамовиц---Стиган, раздел 17.7. Но по прежнему не понял, выражается ли интеграл третьего рода через первый и второй... Может быть, он линейно не выражается?

peregoudov · Сообщение **peregoudov** » 12 дек 2018, 13:43

Так, я проверил формулы из Абрамовиц---Стиган. Оказывается, они работают, а интегал 3-го рода выражается через интегралы 1-го и 2-го родов! Я по-прежнему не понимаю, как это получается, но Maple подтверждает это численно.

Более-менее стандартное определение интеграла 3-го рода

$\Pi(n,k)=\int_0^1\frac{dt}{(1-nt^2)\sqrt{1-t^2}\,\sqrt{1-k^2t^2}}.$

Заодно определение интеграла 1-го рода

$F(z,k)=\int_0^z\frac{dt}{\sqrt{1-t^2}\,\sqrt{1-k^2t^2}},\quad K(k)=F(1,k).$

Аналогично для интеграла 2-го рода.

Справедливы такие формулы.

Для $0<n<k^2$

$\Pi(n,k)=K(k)+\sqrt{\frac n{(1-n)(k^2-n)}}\left[K(k)E\!\left(\frac{\sqrt n}k,k\right)-E(k)F\!\left(\frac{\sqrt n}k,k\right)\right].$

Для $k^2<n<1$

$\Pi(n,k)=K(k)+\sqrt{\frac n{(1-n)(n-k^2)}}\,\Biggl[\frac\pi2-K(k)E\!\left(\sqrt{\frac{1-n}{1-k^2}},\sqrt{1-k^2}\right)-{}$
${}-E(k)F\!\left(\sqrt{\frac{1-n}{1-k^2}},\sqrt{1-k^2}\right) +K(k)F\!\left(\sqrt{\frac{1-n}{1-k^2}},\sqrt{1-k^2}\right)\Biggr].$

Для $n$

вне интервала $0<n<1$

существуют формулы приведения к этому интервалу. Например, для интересного нам случая $n<0$

$\Pi(n,k)=\frac{k^2}{k^2-n}K(k)-\frac{n(1-k^2)}{(1-n)(k^2-n)}\Pi\!\left(\frac{k^2-n}{1-n},k\right).$

peregoudov · Сообщение **peregoudov** » 14 дек 2018, 14:26

Свой правильный ответ у меня есть, но что же делать с ПБМ? Путь к их формуле был долог и труден.

Прежде всего нужно сделать другую подстановку: подгадать так, чтобы одна из точек ветвления, для определенности $t_0=-b/c-ia/c$

, отображалась в бесконечность, а еще одна, для определенности $t=-1$

--- в нуль. Тогда после преобразования под корнем останется многочлен третьей степени, причем один из корней у него будет $z=0$

. Последующей заменой $z\to z^2$ получим стандартный эллиптический интеграл. В общем, замена теперь имеет вид

$t=\frac{t_0z^2+(1-t_0)/2}{z^2-(1-t_0)/2}.$

После такой замены получается интеграл

$I_{1,0}=\frac1b-\frac1{2\pi b}\oint\frac{ab+i(a^2+bc+c^2)-2ibcz^2}{(b+c+ia-2cz^2)\sqrt{1-z^2}\,\sqrt{(b+c)^2+a^2-4bcz^2}}dz.$

Тут уже вырисовываются модуль и параметр

$k=\sqrt{\frac{4bc}{(b+c)^2+a^2}},\quad n=\frac{2c}{b+c+ia},$

сам же интеграл выражается через стандартные

$I_{1,0}=\frac1b-\frac{ik}{\pi\sqrt{bc}}\left(K(k)-\frac{b-c+ia}b\Pi(n,k)\right).$

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

Прежде всего я обнаружил, что для нашего $n$

работает приведенная выше формула

$\Pi(n,k)=K(k)+\frac\pi2\sqrt{\frac n{(1-n)(n-k^2)}}\left(1-\Lambda_0\!\left(\sqrt{\frac{1-n}{1-k^2}},k\right)\right),$

причем множитель перед скобками равен $\dfrac{2i\sqrt{bc}}{k(b-c+ia)}$ , а аргумент $\Lambda_0$ ---

$u=\sqrt{\frac{b+c-ia}{b-c-ia}}.$

Интеграл принимает вид

$I_{1,0}=-\frac{k(a+ic)}{\pi b^{3/2}\sqrt c}K(k)+\frac1b\Lambda_0(u,k).$

Чтобы выделить отсюда вещественную часть явно, нужно сосчитать $\Lambda_0(u,k)+\Lambda_0(u^*,k)$ , что сводится к суммированию неполных интегралов 1-го и 2-го родов с разными амплитудами. Для таких сумм есть формулы

$F(z_1,k)+F(z_2,k)=F(z,k),\quad E(z_1,k)+E(z_2,k)=E(z,k)+k^2z_1z_2z,$

$z=\frac{z_1\sqrt{1-z_2^2}\,\sqrt{1-k^2z_2^2}+z_2\sqrt{1-z_1^2}\,\sqrt{1-k^2z_1^2}}{1-k^2z_1^2z_2^2}.$

Однако в нашем случае они напрямую неприменимы, потому что оказывается $1-k^2z_1^2z_2^2=0$

и получается неопределенность 0/0. Лопиталя ее по $z_2$

, находим формулу для случая $k=1/z_1z_2$

$z=\frac{z_1^2-z_1^2}{2i\sqrt{1-z_1^2}\,\sqrt{1-z_2^2}}.$

В нашем случае результирующая амплитуда

$v=\frac a{\sqrt{(b-c)^2+a^2}},$

а $\Lambda_0(u,k)+\Lambda_0(u^*,k)=\Lambda_0(v,k)+\dfrac{ka}{\pi\sqrt{bc}}K(k)$ . Окончательно

$I_{1,0}=-\frac{ka}{2\pi b^{3/2}\sqrt c}K(k)+\frac1{2b}\Lambda_0(v,k).$

То есть, собственно говоря, в ПБМ лишние квадратные скобки, а так все правильно. Последний член ( $1/b$

) описывает постоянную добавку к полю внутри магнита, а знак $\pm$ во втором члене описывает ветвление при $v=1$

(загиб контура интегрирования за край разреза).

Портал Естественных Наук, версия 1.1