Аннотация:
Построено новое точное решение задачи о тонкой заряженной нити, параллельной границе раздела диэлектриков, один из которых обладает нелинейной восприимчивостью.
Задачи электростатики диэлектриков с нелинейным откликом имеют большое значение как в теоретическом, так и в прикладном аспекте. Такие среды обладают рядом уникальных свойств по сравнению с другими веществами [1], [2]. Для выявления и исследования свойств таких сред разрабатываются новые вычислительные методы, изучаются новые композитные материалы и наноматериалы на основе изоляторов [3]–[5]. Они находят применение в задачах микроэлектроники, интегральной фотоники и др.
Важной составной частью перечисленных прикладных проблем является задача электростатики о взаимодействии тонкой заряженной нити, расположенной над полубесконечной сегнетоэлектрической средой и параллельной границе раздела сред. Математическая постановка задачи включает уравнение относительно электростатического потенциала по обе стороны от границы раздела, условия сопряжения на границе раздела сред и условия на бесконечности.
Точное решение такой задачи неизвестно. Поэтому традиционно ее решают численно или с помощью методов теории возмущений. При численном решении составляют некоторую разностную схему, которая приводит к системе нелинейных алгебраических уравнений [6]. Однако этот подход сталкивается с рядом трудностей.
Во-первых, сеточная реализация условий на бесконечности представляет проблему. Упрощение этих условий и перенос границы раздела на конечное расстояние от нити вносит физическую погрешность в решение, оценить которую затруднительно. Во-вторых, решение указанной нелинейной алгебраической системы сопряжено с большими объемами вычислений. В-третьих, источником поля является заряженная нить, т. е. двумерный точечный заряд. Поэтому решение имеет сингулярность в точке расположения заряда. Такие решения исключительно трудны для численного расчета. Вычитание сингулярности из решения приводит к модификации условий сопряжения на границе раздела, т. е. если источником поля становится не точечный заряд в верхней полуплоскости, а фиктивные поверхностные заряды на границе раздела сред. Такая задача еще труднее для численного расчета.
Методы теории возмущений имеют свои ограничения. Они обеспечивают разумную точность при не слишком больших полях, т. е. если отличие искомого решения от решения линеаризованной задачи достаточно мало. Кроме того, не удается точно удовлетворить условиям сопряжения на границе раздела сред [7]. Это ограничивает физическую точность такого подхода.
В настоящей работе построено точное решение данной задачи в квадратурах, которое сводится к расчету нескольких одномерных интегралов. Это решение справедливо как при малых полях (когда нелинейный отклик среды мал), так и при больших полях (при которых нелинейность становится существенной). Это решение является новым.
2. Постановка задачи
2.1. Физическая постановка
Пусть бесконечно тонкая и бесконечно протяженная заряженная нить с постоянной линейной плотностью заряда $\gamma $ расположена на расстоянии $h$ от плоской границы раздела двух полупространств, заполненных веществом с диэлектрическими проницаемостями $\varepsilon_{1} =1$ и $\varepsilon_{2} =\varepsilon (\mathbf{E})$ соответственно (см. рис. 1).
Среда, заполняющая нижнее полупространство, является сегнетоэлектриком. Для этой среды зависимость вектора поляризации $\mathbf{P}$ от напряженности электрического поля $\mathbf{E}$ является нелинейной. Она имеет следующие качественные особенности. Во-первых, вектор $\mathbf{P}$ сонаправлен вектору $\mathbf{E}$. Во-вторых, при малых полях $E$, много меньших некоторого критического значения $E_\mathrm{cr}$, отклик среды является практически линейным. В-третьих, по мере увеличения поля при $E \sim E_\mathrm{cr}$ начинает проявляться насыщение. При $E \gg E_\mathrm{cr}$ отклик среды перестает зависеть от внешнего поля, и поляризация становится практически постоянной $P \approx P_\mathrm{s}$, где $P_\mathrm{s}$ – поляризация насыщения. Примером такой зависимости $\mathbf{P}(\mathbf{E})$ является следующее выражение:
Возможны и более сложные варианты зависимости $\mathbf{P}(\mathbf{E})$.
Требуется найти распределение потенциала в каждой области в широком диапазоне значений $E$: меньших $E_\mathrm{cr}$, сопоставимых с $E_\mathrm{cr}$ и больших $E_\mathrm{cr}$.
2.2. Система уравнений
Поскольку нить является бесконечно протяженной, решение не зависит от координаты вдоль нити. Достаточно рассмотреть одно выделенное сечение, т. е. задача является двумерной. Введем декартову систему координат, как показано на рис. 1. Тогда границей раздела двух сред является прямая $y=0$, а заряд (выделенное сечение заряженной нити) находится в точке $M_{h} =\{0,h\}$.
Как известно, двумерная задача электростатики описывается двумя уравнениями Максвелла [8]
где $\mathbf{D}=\varepsilon \mathbf{E}$, $\rho$ – плотность свободных зарядов. Отметим, что множитель $2\pi$ обусловлен двумерностью задачи.
Введем скалярный потенциал $\phi$ такой, что $\mathbf{E}= - \nabla \phi$. Подставим последнее выражение в (2), учитывая, что в верхнем полупространстве присутствует заряженная нить, а в нижнем полупространстве заряды отсутствуют. Получим уравнения электростатики
Условие (5) вытекает из непрерывности потенциала в любой точке пространства. Условие (6) означает, что единственной причиной разрыва нормальной составляющей вектора электрической индукции $\mathbf{D}$ на границе раздела двух сред является наличие на ней поверхностных свободных зарядов. Перепишем уравнение (6) для потенциалов, учитывая, что свободных зарядов на границе нет, а нормаль сонаправлена с осью $y$, которая изображена на рис. 1. Получим
Поскольку задача рассматривается в неограниченной области, необходимо поставить условия на бесконечности. По мере удаления точки наблюдения от заряда напряженность поля ослабевает. Вклад нелинейности нижней среды уменьшается, и диэлектрическая проницаемость стремится к постоянной величине, которую мы укажем далее (см. п. 3.2). Поэтому условие на бесконечности имеет вид
при $\rho =\sqrt{x^{2} +(y-h)^{2} } \to \infty $. Здесь $\phi ^{0} $ – решение линейной задачи, в которой диэлектрическая проницаемость нижней среды является постоянной величиной.
Отметим, что такие условия на бесконечности являются нетривиальными. Они не сводятся к классическим условиям регулярности (ограниченность либо стремление к нулю). Условия (8) ранее в литературе не предлагались.
Таким образом, окончательная постановка задачи состоит из уравнений (3)–(5), (7), (8).
3. Нижняя полуплоскость
Основную трудность представляет построение решения в нижней полуплоскости. Сначала мы рассмотрим предельные случай, когда это решение удается построить в явном виде. Затем обратимся к общему случаю.
3.1. “Большие” поля
Под “большими” полями мы подразумеваем предельный случай, когда $E\gg E_\mathrm{cr}$. Тогда экспонентой в формуле (1) можно пренебречь, а поляризация будет приближенно равна поляризации насыщения
Решение уравнения (10) должно быть центрально-симметричным относительно точки $M_{h} =\{0,h\}$, поскольку единственным источником поля является заряд в указанной точке. Введем полярную систему координат с центром в точке нахождения заряда, как показано на рис. 1:
В этих координатах $\mathbf{E} = E \mathbf{e}_{\rho}$, $\mathbf{P}=P_\mathrm{s} \mathbf{e}_{\rho}$. В двумерном случае $\mathbf{D}= \mathbf{E} + 2\pi \mathbf{P}$, поэтому
где $A$ – константа интегрирования, не зависящая от координаты $\rho$. Мы найдем ее позже (см. п. 3.4). Отметим, что уравнение для индукции и его решение не зависят от свойств среды.
Подставляя (14) в (12), получим напряженность поля $E$ в нижней полуплоскости в пределе “больших” полей:
Другим предельным случаем являются “малые” поля. Под ними будем понимать случай, когда $E\ll E_\mathrm{cr}$. В этом случае экспоненту в (1) можно разложить в ряд Тейлора и ограничиться линейным членом. Тогда вектор поляризации имеет вид
Отсюда видно, что в случае малых полей диэлектрическая проницаемость нижней среды является константой и не зависит от поля.
Решение этой линейной задачи строится методом изображений и хорошо известно. Пусть $M=\{x,y\}$ – точка наблюдения, $M_{-h} =\{0,-h\}$ – изображение точки $M_{h} $ относительно прямой $y=0$. Потенциал электрического поля равен
Величины $\varepsilon_{1} $ и $\varepsilon_{2} $ суть диэлектрические проницаемости верхней и нижней сред соответственно. В нашем случае $\varepsilon_{1} =1$, а $\varepsilon_{2}$ определяется формулой (18). Напряженность поля имеет следующий вид:
Теперь обсудим общий случай, не сводящийся к рассмотренным выше. Пусть характерная величина поля $E$ сопоставима с $E_\mathrm{cr} $. Построим решение аналогично п. 3.1. Однако надо учесть, что вектор поляризации описывается формулой (1) без каких-либо пренебрежений. В результате приходим к следующему алгебраическому уравнению относительно $E$:
Уравнение (23) является трансцендентным и хорошо решается численно, например, с помощью метода Ньютона. Решив (23) и найдя значения поля, мы можем найти потенциал, как при получении (16).
Задачи с более сложными зависимостями $P(E)$ рассматриваются аналогично. При этом случаи сильных и промежуточных полей полностью повторяют п. 3.1 и 3.3, а для слабых полей может отличаться конкретное выражение для $\varepsilon_2$.
Таким образом, мы построили решение задачи в нижней полуплоскости в квадратурах. Осталось только найти константы интегрирования $A$, $B$, $C$.
3.4. Константы интегрирования
Найдем константы интегрирования из условий на бесконечности. При достаточно большом $\rho \to \infty$ решение для “промежуточных” полей должно переходить в решение для “слабых” полей. В частности, потенциал должен стремиться к выражениям (19), (20), а напряженность – к (21), (22).
Константу $C$ найдем, сравнивая напряженности “промежуточного” и “слабого” полей. В формуле (23) будем считать поле $E$ достаточно малым. Тогда поляризация определяется формулой (17). Подставим (17) в (23) и решим полученное уравнение относительно $E$. Получим
Теперь определим константу $B$, сравнивая нормировку потенциалов “промежуточного” и “малого” полей. В ходе расчета потенциала “промежуточного” поля в нижней среде мы вычисляем интеграл с переменным верхним пределом от напряженности по переменной $\rho$, причем нижним пределом является $\rho_\mathrm{min} =h$, верхним – текущее $\rho_{i}$:
На рис. 2 приведен график зависимости модуля разности решений для “промежуточного” поля и для “малого” поля (сплошная линия) и производной этой разности (штриховая линия) от координаты $\rho$ в двойном логарифмическом масштабе. Из графика видно, что разность решений выходит на константу, которая и является константой нормировки. Добавляя эту константу к (26), получим окончательное решение для “промежуточных” полей в нижней среде.
Найдем теперь константу $A$ в решении для “больших” полей. Очевидно, что при малых $\rho$ решение уравнения (23) должно переходить в решение (15). В обоих случаях преобладающим является слагаемое $\sim 1/\rho$. Отсюда заключаем, что $A=C$.
3.5. Пример расчета
Пусть $\varepsilon_{2} =2$, $E_\mathrm{cr} =0.05$, $\gamma =1.5$, $h = 1$. На рис. 3 показаны потенциалы для случаев “малых”, “промежуточных” и “больших” полей в нижней полуплоскости в зависимости от $\rho$. Видно, что при $\rho \to h$ потенциал “промежуточного” поля переходит в потенциал для случая “больших” полей. Это соответствует теоретическим ожиданиям. При значениях $\rho <3$ потенциалы “большого” и “промежуточного” полей практически совпадают. В этой области нижней среды имеет место насыщение.
За пределами области насыщения с увеличением $\rho$ выражение для “больших” полей довольно быстро теряет точность. Начиная с некоторого $\rho$ (в данном примере $\rho \approx 7$), потенциал “большого” поля начинает возрастать, т. е. применение данного предельного случая становится неправомерным.
Видно также, что при достаточно больших $\rho $ потенциал “промежуточного” поля переходит в решение для “малого” поля. Это соответствует теории. Для тех $\rho$, при которых потенциалы “малого” и “промежуточного” полей визуально совпадают, нижнюю среду можно с хорошей точностью считать линейной.
При малых $\rho$ решение для “малого” поля качественно воспроизводит решение для “промежуточного” поля, однако с уменьшением $\rho$ расхождение увеличивается. Это объясняется тем, что потенциалы “большого” и “малого” полей имеют разные коэффициенты перед логарифмическим слагаемым.
Таким образом, проведенный расчет демонстрирует правильность предельных переходов потенциала “промежуточного” поля при больших и малых $\rho$. В частности, это показывает правильность выбора постоянных интегрирования.
4. Верхняя полуплоскость
4.1. Третья формула Грина
Решив задачу в нижней полуплоскости, мы можем вычислить потенциал и его производную на границе раздела. Для нахождения потенциала в верхней полуплоскости можно воспользоваться третьей формулой Грина. В двумерном случае она имеет вид [9]
Здесь $D$ – расчетная область (в нашем случае $y > 0$), $L$ – ее граница ($y = 0$), $P = \{x ,0\}$ – точка на границе, $M_{0} =\{x_0, y_0\}$ – точка наблюдения, $n_{P} $ – нормаль к границе,
В первом слагаемом в (29) учтем, что нормаль к границе раздела сонаправлена с осью $y$, а граница раздела является прямой $y=0$. Тогда это слагаемое принимает вид
Отметим, что решение в верхней полуплоскости не является центрально-симметричным. Это непосредственно видно даже в линейном случае, поскольку центральная симметрия нарушается из-за наличия заряда-изображения в точке $M_{-h} $. При этом сохраняется осевая симметрия относительно прямой $x=0$.
Интегрирование в (36) ведется по всей числовой прямой. Для удобства введем замену переменных $x = \xi (1 - \xi^{2})^{-m}$, после которой интеграл преобразуется к следующему виду:
В дальнейших расчетах было выбрано значение $m=1.1$. Замена переменных (37) сводит интегрирование по всей прямой к интегрированию по отрезку $[-1,1]$, на котором уже можно ввести равномерную сетку и использовать классические квадратурные формулы.
4.2. Пример расчета
Проверим справедливость решения (36). Для этого рассмотрим случай “малых” полей и вычислим потенциал верхней среды двумя способами. С одной стороны, этот потенциал выражается формулой (20). С другой стороны, мы можем вычислить его по формуле (36), подставляя в нее потенциал нижней среды (19). Полученные решения должны совпасть с точностью до погрешности численного интегрирования.
Для $\varepsilon_{2} =2$ она приведена на рис. 4 в виде изолиний с фиксированным шагом. Цифры на изолиниях – соответствующие значения представляемых величин. Видно, что относительная разность практически во всей расчетной области не превышает 0.1%. Это соответствует точности квадратурной формулы на выбранной сетке. Таким образом, данный расчет показывает справедливость подхода к расчету поля в верхней среде, основанного на третьей формуле Грина.
Напомним алгоритм решения исходной нелинейной задачи. Уравнение (23) численно решается методом Ньютона, в результате чего получается массив сеточных значений поля в нижней среде при различных сеточных значениях $\rho$. С помощью (27) вычисляются значения потенциала нижней среды при тех же значениях $\rho$. Затем значения потенциала на границе подставляются в интеграл (36), который рассчитывается численно с заменой переменных (37). После данной процедуры мы получаем потенциал поля в различных точках $(x,y)$ как выше, так и ниже границы раздела.
В качестве примера проведем расчет при $h=1$, $\gamma =10$, $E_\mathrm{cr} =0.05$ и проиллюстрируем полученное решение. На рис. 5 показано распределение потенциала в пространстве в виде изолиний. В точке расположения заряда потенциал имеет сингулярность. Вблизи этой точки значения $\phi$ быстро нарастают. Видно также, что на границе раздела поле преломляется. Потенциал непрерывен, а его нормальная производная испытывает разрыв.
Проиллюстрируем роль нелинейности. Вычислим модуль электростатического поля и сравним с критическим полем $E_\mathrm{cr}$. На рис. 6 показано распределение поля $E$ в пространстве в виде изолиний. Видно, что в нижней полуплоскости на достаточно большом удалении от начала координат имеем $E \ll E_\mathrm{cr}$. На рис. 6 имеется область $E/E_\mathrm{cr} < 0.5$, в этой области нелинейные поправки малы. При $\rho \sim 25$ величина $E$ становится сопоставимой с $E_\mathrm{cr}$, и начинают проявляться нелинейные эффекты. Наконец, при $\rho < 2.5$ имеем $E \gg E_\mathrm{cr}$. В этой области отклик среды является сильно нелинейным, а поляризация близка к $P_\mathrm{s}$. В верхней полуплоскости поле $E$ имеет сингулярность. На границе раздела поле преломляется и его модуль испытывает разрыв.
Таким образом, при выбранных параметрах решение не сводится к какому-то одному предельному случаю, а охватывает случаи $E \ll E_\mathrm{cr}$, $E \sim E_\mathrm{cr}$, $E \gg E_\mathrm{cr}$. Поэтому данный пример достаточно представителен.
Оценим количественно влияние нелинейности на решение. Для этого сравним решение $\phi$ исходной нелинейной задачи и решение $\phi_\mathrm{lin}$ линеаризованной задачи (см. п. 3.2) при тех же значениях $h$, $\gamma$. На рис. 7 показана относительная разность этих решений
Нормировка в формуле (39) выбрана из следующих соображений. Потенциалы исходной и линеаризованной задач меняют знак в пределах расчетной области. Поэтому нормировать ни на тот, ни на другой нецелесообразно, поскольку при делении может получиться сколь угодно большая величина. Однако геометрическое место точек смены знака для этих потенциалов различно. Поэтому сумма модулей $|\phi_\mathrm{lin}| + |\phi|$ не обращается в нуль.
Из рис. 7 видно, что пределом решения при радиусе, стремящемся к бесконечности, являются выражения (19), (20). Это подтверждает, что при достаточном удалении от заряда вклад нелинейности ослабевает. По мере приближения к заряду относительная разность (39) увеличивается и достигает $\sim 1.6$. Тем самым при выбранных значениях параметров вклад нелинейности нижней среды оказывается существенным.
В непосредственной близости от заряда расхождение решений быстро уменьшается. Изолинии величины (39) имеют вид почти концентрических окружностей. Наибольшая из них соответствует уровню $1.6$, следующая – $0.8$ и т. д. Это объясняется тем, что решения $\phi$ и $\phi_\mathrm{lin}$ имеют одинаковую асимптотику вблизи сингулярности. Это непосредственно видно из формул (36) и (22).
Конфликт интересов
Авторы заявляют, что у них нет конфликта интересов.
Список литературы
1.
T. Wei, K. Liu, P. Fan et al., “Novel NaNbO$_3$–Sr$_{0.7}$Bi$_{0.2}$TiO$_3$ lead-free dielectric ceramics with excellent energy storage properties”, Ceram. Int., 47:3 (2001), 3713–3719
2.
J. Shi, A. H. Akbarzadeh, “Hierarchical cellular ferroelectric metamaterials: A design motif to enhance multifunctional figures of merit”, Comp. Struct., 250 (2020), 112395, 20 pp.
3.
C. C. Wang, G. Pilania, S. A. Boggs, S. Kumar, C. Breneman, R. Ramprasad, “Computational strategies for polymer dielectrics design”, Polymer, 55:4 (2014), 979–988
4.
J. Liu, M. Zhang, L. Guan et al., “Preparation of BT/GNP/PS/PVDF composites with controllable phase structure and dielectric properties”, Polymer Testing, 100 (2021), 107236, 10 pp.
5.
L. Chen, L. Yan, Y. Guo et al., “Chemically functionalized multi-walled CNTs induced phase behaviors of poly (Vinylidene Fluoride) nanocomposites and its dielectric properties”, Synthetic Metals, 269 (2020), 116555, 10 pp.
6.
Н. Н. Калиткин, П. В. Корякин, Численные методы, т. 2, Методы математической физики, Академия, М., 2013
7.
М. Л. Акимов, М. А. Пятаков, О. П. Поляков, П. А. Поляков, “Электростатические поля поляризационных зарядов в диэлектрике с нелинейной восприимчивостью”, Изв. РАН. Сер. физ., 86:9 (2022), 1251–1255
8.
В. И. Денисов, Введение в электродинамику материальных сред, Изд-во Моск. ун-та, М., 1989
9.
А. Н. Тихонов, А. А. Самарский, Уравнения математической физики, Изд-во Моск. ун-та, М., 1999
Образец цитирования:
А. А. Белов, М. А. Тинтул, П. А. Поляков, “Точное решение задачи о взаимодействии точечного заряда и диэлектрика с нелинейной восприимчивостью”, ТМФ, 220:3 (2024), 578–590; Theoret. and Math. Phys., 220:3 (2024), 1556–1566