ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
2012
Серия: Физика
Вып. 3 (21)
УДК 537.2
Взаимодействие заряженных проводников на близких расстояниях
Е. Л. Тарунин
Пермский государственный национальный исследовательский университет, 614990, Пермь, ул. Букирева, 15
Методом сеток решена задача о взаимодействии двух заряженных проводников в форме цилиндров (шайб). Распределение потенциала вне проводников находилось из решения соответствующей задачи Дирихле. При программировании использовалась версия языка Паскаль АВС, позволяющая использовать массивы большого размера. Подробно рассмотрен случай взаимодействия одноименных зарядов. Как и в [1] показано, что при малых расстояниях между телами сила взаимодействия существенно отличается от закона Кулона для точечных зарядов.
Ключевые слова: закон Кулона, метод сеток, задача Дирихле.
1. Введение
При сближении заряженных проводников происходит перераспределение зарядов на поверхности проводников. Это перераспределение заряда меняет в свою очередь и силу, действующую на проводники. В работах [1-3] показано, что при сближении шаров отклонение от закона Кулона весьма существенно. В данной работе рассмотрено взаимодействие заряженных проводящих тел в форме цилиндров (шайб). Аналогичный подход использован в работе [4] при вычислении сил, действующих на заряженные тела в поле плоского конденсатора. При программировании использована версия языка Паскаль АВС, позволяющая использовать массивы большого размера.
2. Постановка задачи и метод решения
Геометрия цилиндрической области изображена на рис.1.
R3 V 1 1 1 1 1
R2 IV 1 1 1
R1 I II 1 III 1
0 x x x x
Рис.1. Геометрия области
x
Круглые шайбы расположены на оси, проходящей через их центры. Здесь Я1 - радиус первой шайбы, - радиус второй шайбы (Я2 — Я1). Толщины шайб равны соответственно Ь1=х2-х1, Ь2=х4-х3. Римскими цифрами обозначены области интегрирования. При равенстве радиусов шайб область IV отсутствует.
Использованы основные положения электростатики [5]: свободные заряды на проводящих телах располагаются только на его поверхности, поверхностная плотность заряда на проводнике пропорциональна нормальной составляющей индукции электрического поля
ст = Dn = ssoЕп .
(1.1)
Здесь s - диэлектрическая постоянная среды,
s0 - const ( s
Напряженность электрического поля связана с потенциалом поля соотношением
¡8.85 -10~12Кл2/(Н • м2).
E = -grad (<р).
(1.2)
Полный заряд на проводнике вычисляется по формуле
2 = | Б^Б = 8 | Е,¿Б = ~ее0 ¡^ dS. (1.3)
Задача решается в безразмерных переменных. Полагается, что 8 = 1, 80 = 1. В качестве единицы расстояния выбран радиус первой шайбы Я1.
Расчеты выполняются в предположении осевой симметрии в области цилиндрических координат (х, г). Уравнение Лапласа, записанное в цилиндри-
© Тарунин Е. Л., 2012
r
ческих координатах для потенциала, имеет вид
. 1 д др д р
Ар =--(г —) = 0.
г дг дг дх
(1.4)
Перейдем к описанию параметров задачи и метода сеток. Для расчетов используется квадратная сетка х ='• к, гк = к • к. Параметрами задачи являются: радиусы шайб R1, Я2, толщины шайб Ь1, Ь2, расстояние между ближайшими торцами шайб Ь и заряды на шайбах Q1, Q2. Параметрами метода являются: к - шаг пространственной сетки, размеры области, на которой задается распределение потенциала, точность решения уравнения Лапласа итерационным методом. Разностная аппроксимация уравнения Лапласа (1.4) на 5-ти точечном шаблоне [6,8] имела вид
У[1,к] = (С1[к] • V[г,к +1] + С2[к] • V[г,к -1] + + V [г +1, к] + У[1 -1, к])/4.
(1.5)
в которой г}- (] = 1,2) - расстояния до центра соответствующей шайбы с потенциалом Vsj были увеличены так, чтобы в точках, близких к шайбам, значения потенциала были меньше значения потенциала на шайбе.
Итерации прекращались при выполнении неравенства (максимум определялся по всем внутренним узлам области):
тах
К[г,к] - V(т)[г,к]| < е = 36к2е0
е0 = 10-
.(1.9)
При записи (1.5) заменено обозначение потенциала (р ^ V). Уравнение приведено к виду, удобному для использования итерационного метода. В целях сокращения затрат машинного времени переменные коэффициенты в (1.5) вычислялись до начала итераций С1[к]=1+0.5/к, С2[к]=1-0.5/к.
Для узлов на оси области для итерирования использовалась формула
V[г,0] = (V[г +1,0] + V[г -1,0] + 4 • V[г,1])/6, (1.6)
которая получается при раскрытии неопределен-
г /1 дрч
ности типа ноль на ноль предела пт(--—) при
г дг
г ^ 0 по правилу Лопиталя. Погрешность использованной аппроксимации уравнения Лапласа О( к2).
По методу Зейделя, примененного к формулам (1.5), (1.6), определялись предварительные значения V[г,к]. Окончательные значения на (т +1) итерации вычислялись с использованием идеи поточечной верхней релаксации [7]:
V(m+1)[/, к] = V(m)[/, к] + ю- (Гф, к] - V(m)[/, к]). (1.7)
Оптимальный параметр верхней релаксации ф определялся по формуле [7], зависящей от числа интервалов по осям прямоугольной области интегрирования х5,Я3. Итерации выполнялись с левого нижнего узла области (1=1, к=1) направо и снизу вверх.
Начальное значение потенциала во внутренних узлах области за пределами шайб равнялось нулю или вычислялось по формуле
V (°) = ™+ г1 г2
(1.8)
Часть расчетов выполнялась с уточнением сеточных значений потенциала, полученных на пятиточечном шаблоне. Уточнение выполнялось с помощью дополнительных итераций по методу Зейделя для разностных уравнений, полученных на девятиточечном шаблоне. Так как введение дополнительных итераций практически не давало уточнения интегральных характеристик, от них пришлось отказаться.
Задание числа интервалов на радиусе первой шайбы N0 (10 < N0 < 30) определяло шаг пространственной сетки к = 1/ N0 и целочисленные значения индексов для всех граничных значений (напомним, что все расстояния являются относительными по отношению к радиусу первой шайбы Я). Внешние границы области (х = 0, х = х5, г = Я3) отдалены от шайб на расстояние, пропорциональное ¡л^ Я (л - множитель, параметр метода). Внешний радиус области Я3 вычислялся по формуле Я3 = ¡Я1. Левая и правая границы области отстояли от ближайшего торца на расстоянии х = л • кх • Я . Параметр кх так же, как и л , естественно влиял на результаты расчетов. Определение "оптимальных" значений этих параметров выполнялось с помощью вычислительных экспериментов (ВЭ). Основные расчеты выполнялись при значениях л = 40, кх = 0.6 . Чтобы учесть особенность угловых точек, часть расчетов выполнялась для тел с округленными углами (радиус округления равнялся шагу сетки).
Опишем алгоритм определения потенциалов V, ^ на поверхностях шайб при задании на них зарядов. Такую задачу называют "задачей Робена" [8]. В соответствии с электростатикой предполагалось, что заряды на проводниках определяются их линейной зависимостью от их потенциалов:
J
ъ =Е С ^. (1.10)
j=1
Здесь Сг,j - емкостные коэффициенты; Сгг -коэффициент электроемкости г - го тела (Сй > 0); Сг j - коэффициент индукции между г -м и j -м телами, когда потенциал j -го тела равен Vj, а ос-
тальные тела имеют потенциал равный нулю
Су = С^ < 0).
Для нахождения емкостных коэффициентов для двух тел (J = 2) решаются две вспомогательные задачи. Верхними индексами будем отмечать номер задачи. В первой задаче потенциалы тел задавались равными V(1) = 1, К2(1) =0, а во второй задаче V(2) =0, К2(2) =1. После решения соответствующих задач Дирихле и нахождения зарядов q(m')(i = 1,2; m = 1,2) по формулам, вытекающим из (1.10), определяем коэффициенты системы
О = ф(1) • ф«2) - ф» • q1(2),
О = Ф? • q\ - 4т1 • q2,
О2 = Ф0) • Ф2 - Ф12) •
(1.11)
(1.13)
мах (1.13) использовались значения (гу- • /у).
Элемент силы, действующий на элемент поверхности х = х, вычислялся по формуле [5]:
/ =-0.5 • Е2т • dS у.
(1.14)
и по правилу Крамера определяем значения потенциалов.
Решив третью задачу с вычисленными потенциалами V, У2, определяем распределение зарядов на телах и вычисляем силу взаимодействия между ними.
Опишем формулы вычисления распределения зарядов, используемые во всех трех задачах. С учетом принятых допущений (е = е0 = 1) поверхностная плотность заряда стп = Еп = . Полный
дп
заряд представляем в виде суммы трех слагаемых (аналогичные формулы используются для второго заряда) для трех поверхностей при х = х1, х = х2 и г = Я1:
д1=0х1+0х2+0т1. (1.12)
Производная по нормали вычислялась по разностным отношениям - со вторым порядком аппроксимации 0(к2). В одномерном случае аппроксимация со вторым порядком имеет вид
2 V -0.5 V -1.5•V0^dV_h2 dV к dx 3 скъ
Интеграл по поверхности вычислялся по формуле Симпсона с остаточным членом :
| / (х^х = к 2 (/у-1 + 4 /у + /у +1),
Я =- — Г (ш) Я = 9 1 .
Суммирование в (1.15) идет только для нечетных индексов у . Поэтому число интервалов на интегрируемом участке должно быть четным. При вычислении поверхностных интегралов на торцевых поверхностях элемент поверхности равнялся dS = 2жdr . Поэтому вместо значений /■ в сум-
Знак "минус" в этой формуле обусловлен тем, что внешняя нормаль к поверхности х = х противоположна положительному направлению оси х . Аналогичная формула использовалась для вычисления силы Г12, действующей на торец шайбы при х = х . В итоге сила, действующая на первую шайбу, вычислялась как сумма Г = Гп + Г12. Заметим, что сила, действующая на боковую поверхность шайбы (г = Я1), не дает вклада в компоненту силы в направлении оси х .
Основное время при расчетах тратится на итерационный метод решения задачи Дирихле. Объем вычислений зависит от числа узлов области и от требуемой точности. Число требуемых итераций в
каждой задаче при е1 = 10-14 обычно не превосходило т = 8 • Ыг 3 .
Использовалось два варианта граничных условий. В первом варианте значение потенциала на внешней границе задавалось равным нулю (этот вариант будет обозначаться как Vg = 0), а во втором варианте, обозначаемом как Vg = 1, значения потенциала вычислялись по формулам (1.8), которые использовались и при вычислении начальных значений.
3. Характеристики решения
Основной характеристикой решения задачи является безразмерный коэффициент, определяемый через отношение сил
к = г,/^ . (2.1)
Здесь Гу - сила, действующая на первый заряд и вычисленная из решения основной (третьей) задачи, а Г0 - сила, вычисленная по формул4)Кулона в предположении, что заряды сосредоточены в центре шайб:
/т б1* 62 61* 62
Г0 =-=-. (2.2)
4ж22 4Л(1 + 0.5*( Ц + Ь2 ))2
Коэффициент к показывает относительное отклонение от закона Кулона. Из постановки задачи следует, что этот коэффициент зависит от большого числа параметров. В этой статье рассматривается вариант с равными зарядами, которые располагаются на шайбах с равными значениями радиусов и толщин: Р1=Р2, Я1=Я2= Ь1=Ь2. В этом случае основным параметром, определяющим значение этого коэффициента, является расстояние между ближайшими торцами шайб Ь.
Величина силы взаимодействия ^ определялась интегрированием по формулам Симпсона для подынтегральных выражений типа (1.16). В части расчетов эта сила вычислялась через производные от потенциальных энергий
F = -
~дЬ
F = --öW
dL
(2.3)
Потенциальные энергии системы зарядов вычислялись по формулам [5]
W = 0.5(Vs1Q1 + Vs2 Q2), W = 0.5JE2dV .
(2.4)
кз =
max E„
dQx = Qx1 - Qx2,
(2.6)
Точность решения задачи Дирихле проверялась по выполнимости теоремы Гаусса [5]:
\EndS = q1 + q2 .
(2.7)
Интеграл по замкнутой поверхности в (2.7) вычислялся методом Симпсона с порядком аппроксимации компоненты напряженности О( к2 ) . В интеграл входили соответствующие суммы: ^ (от двух торцов при х = к и х = х5 - к) и (поток через поверхность г = Я3 - к). Отклонение от выполнимости теоремы Гаусса в процентах вычислялось по формуле
Использование формул (2.3) требует счета задач с близкими значениями Ь , так чтобы можно было вычислить производные с помощью конечных разностей с малой погрешностью. Заметим, что все три подхода к вычислению силы взаимодействия дали близкие значения (при отсутствии погрешностей эти значения должны быть равными).
Плотность распределения зарядов на поверхности проводников пропорциональна нормальной компоненте напряженности поля. В расчетах вычислялась зависимость Еп (я) =ст(^) и строился соответствующий график. Здесь 8 - расстояние от угловой точки первой шайбы ( х ,0 ) до соответствующей точки на шайбе. Зависимость Еп (5) = ст(^)
для Ь = 2 представлена на рис. 3 . Как видно, локальные максимумы напряженности поля достигаются в углах области. Интегральной характеристикой обсуждаемой зависимости является коэффициент
DG = 100 • (
S1 + S2 - 1) . qi+q2|
(2.8)
При щ + д2\ =0 (это соответствует зарядам
разных знаков, равных по модулю) эта формулам не использовалась. Заметим, что в этом случае симметрия задачи не дает отклонения от теоремы Гаусса.
Расчеты показали, что при нулевых значениях потенциала на внешней границе (вариант Vg = 0)
значения коэффициента к были большие, чем в варианте Vg ф 0 . Различие соответствующих значений давало косвенную информацию о точности расчетов.
4. Результаты расчетов
Рассмотрим результаты расчетов для фиксированного набора параметров
qi = q2 =1 > Ri = L =Ri= l2 =1.
(3.1)
• »г (2.5)
min En
При уменьшении расстояния между шайбами величина этого коэффициента возрастает (к3 « 3.74 при L = 10 и к3 и 20.5 при L = 1). При использовании варианта расчета с округлением углов зависимость En (s) = ct(s) сглаживается вблизи углов, но величина коэффициента къ меняется несущественно.
Интегральными характеристиками распределения зарядов являются доли зарядов на трех поверхностях: Qx1 при x = x1, Qx2 при x = x2 и QR1 при r = R1 (при L =10 , Vg = 0, например, Qx1 и 22.07% , Qx2 и 21.32% , QR1 и 56.62%).
Вычислялась также разница долей зарядов на торцевых поверхностях:
которая характеризует взаимное влияние заряженных проводников друг на друга. При сближении шайб эта разница возрастает.
Этот набор параметров описывает взаимодействие шайб с одинаковыми размерами ( Я = Я ,
Ь1 = Ь2 ), на которых задан одинаковый заряд. Расстояние между ближайшими торцами шайб изменялось от 10 до 0.25. Параметры л и кх, определяющие расстояния до внешней границы, менялись лишь в предварительных расчетах (обычно л = 40 , кх = 0.6 ).
В случае Vg = 1 значения потенциала на внешней границе области вычислялись по формуле (1.8). Максимальное значение потенциала на границе для типичных значений л, кх достигалось в координате (0,0). Это значение потенциала при Ь = 10 было в 14.5 раз меньше значения потенциала на шайбе. При продвижении по радиусу вдоль левой границы значение потенциала убывало до наименьшего значения, достигаемого в координате (0, Я3). При продвижении вдоль верхней границы (г = Я3) потенциал возрастал до локального максимума, достигаемого в координате (х2 + Ь /2, Я3 ),
а затем убывал в соответствии с симметрией для равных зарядов. Локальный максимум был в 1.456 раз меньше максимального значения. При значении кх и 0.9016 (Ь = 10) локальный максимум оказывается равным глобальному максимуму.
Обсудим сначала процедуру и результаты расчета для набора параметров (3.1) при Ь = 10, N0 = 20, Vg = 0 , е = 9 • 10 . При выбранном значении N 0 = 20 три этапа задачи решаются на пространственной сетке размером 1200*800. Для слежения за состоянием счета на экран выводились: этап решаемой задачи, номер итерационного приближения т с шагом Дт =100 и текущее значение невязки пву. При выбранном шаге для вывода Дт =100 на экране появлялись новые значения примерно через 10 секунд. При первом выводе (т = 100) значение невязки равнялось пеу=0.0408. Затем невязка монотонно убывала, но на интервале т от 1100 до 1700 ее значение возросло от
1.11 •Ю-5 до 4.03 •Ю-5. Эффект временного возрастания невязки типичен для метода ПВР на сетке с большим числом узлов. При т1 =4646 была достигнута требуемая точность решения первого этапа. Полное число итераций на всех этапах 8ш=13 621.
После окончания счета выведены характеристики решения:
ОО = 6.19%, к = 1 1794,
Vs1 = Vs2 = Ж = 0.0906446,
Ж1 = 0.0961106, к3 = 3.623, к5 = 1.9625%, (3.2)
6х1 = 22.066%, 6г1 =56.622 % , dQ = 3.419 .
По значению коэффициента к1 видно, что отклонение от закона Кулона при Ь = 10 незначительно (уточнение дает еще меньшее значение 1.12175). Отличие Ж1 от Ж (в данном расчете около 6%) свидетельствует о погрешности вычисления этих значений; их различие уменьшается при уменьшении шага сетки. Любопытно, что наибольшая доля заряда (более половины) приходится на боковую поверхность при г = Я1. Разница долей зарядов на торцевых поверхностях dQ позволяет понять причину уменьшения сил отталкивания заряженных проводников. При сближении проводников эта величина увеличивается.
Распределение потенциала показано на рис. 2 для части левой половины области. На осях указаны номера узлов. Вторая половина области восстанавливается из соображений симметрии относительно середины между шайбами х = х2 + 0.5Ь . Заметно быстрое убывание потенциала по мере удаления от шайбы.
Рис. 2. Изопотенциальные линии при Ь
Представление о распределении нормальной компоненты напряженности поля Еп (я) на шайбе при Ь = 2 дает рис. 3. Заметим, что в силу связи (1.1) такую же зависимость имеет плотность заряда. Здесь я - расстояние (точнее, номер узла) от левой точки первой шайбы ( х ,0 ) до соответствующей точки на шайбе при продвижении по часовой стрелке.
=10 (Vsl =0.0906446)
Как видно, максимум ст = достигается вблизи углов. Отношение максимального значения max к минимальному min дается коэффициентом къ = 3.623 . Заметно, что напряженность поля в левом углу значительно больше напряженности поля в правом углу. Для расстояния между шайбами L = 10 этот эффект слабее.
150
100-
0
1 т I т I
20 40 60
Рис. 3. Распределение нормальной компоненты напряженности на первой шайбе при Ь = 2
Из теории [9] известно, что на торцевой поверхности заряженного тонкого диска радиуса г = а плотность заряда описывается формулой
а(г) /Vа2 - г 2ла
(3.3)
Практическая зависимость а (г) близка к (3.3), но у края диска устанавливается конечное, хотя и большое, значение, зависящее от особенностей диска. В нашем случае (хотя это не ситуация одного и не тонкого диска) следствие зависимости (3.3) в виде
а (г) = а (0)/7Г - г2 подтверждается примерно до половины радиуса шайбы (значение а(0.9, х) меньше предсказания примерно на 25%).
Обсудим качественное изменение характеристик решения задачи при уменьшенном расстоянии между шайбами Ь = 2 (остальные параметры задачи и метода не изменены). Интегральные характеристики в этом случае таковы:
ПО = 6.221%, к = 1.4617, Vs1 = Vs2 = Ж = 0.108900, Ж1 = 0.115521, к3 = 6.502, к5 = 1.6773%, (3.4) 0>х\ = 23.034%, 0г1 =58.672 % , йЯ = 6.740 .
пространственной сетки к = 1/ N0 . Следующие серии расчетов выясняли зависимость интегральных характеристик, и в первую очередь отношения сил к , от шага сетки.
В табл. 1 даны значения коэффициента к при Ь = 10 , полученные для разных значений N0 = 1/ к при двух вариантах граничных условий Vg = 0 и Vg ф 0. Значения коэффициента, полученные при граничных условиях Vg = 0, обозначены как к10, а значения коэффициента, полученные при граничных условиях Vg ф 0, обозначены как к/.
Таблица 1. Зависимость коэффициента к^0) при Ь = 10
N 0 10 12 16 20
1.19162 1.190452 1.185308 1.179419
к1 1.180962 1.181258 1.178063 1.173441
Дкх 0.010658 0.008194 0.007245 0.005978
продолжение табл. 1
Эти характеристики получены на сетке размером 1040*800 (N0 = 20). По сравнению с (3.2) существенно увеличились: кх, Vs1, Ж1, к3,
ёЯ.
Перейдем к рассмотрению погрешности получаемых решений. Основным параметром метода служило число интервалов на радиусе первой шайбы N0 , определяющее шаг квадратной
N 0 22 24 Нш к
к? 1.176558 1.17381 1.12175
к1 1.171059 1.16972 1.12175
дк 0.005499 0.00509 0.000001
Как видно, при N0 > 12 значения коэффициентов монотонно убывают с ростом N0 . Так же ведет себя и положительная разность значений Дк]= к" - к\. Для относительной разницы получена оценка 5кг и 12 • к %; из нее следует, что при к = 1 / 24 относительное отличие к10 от к11
2
менее 0.5%.
В последнем столбце таблицы указаны предельные значения коэффициентов, полученные с помощью аппроксимации полиномом второй степени, построенном по трем наименьшим значениям шага сетки Р2 (к) = Ишк + с • к + с2 • к2. Соответствующие коэффициенты в полиномах (к < 0.05) таковы:
С0 = 1.72965, с0 = -11.52624,
с\ = 1.594536, с\ = -11.2147. (3.5)
Предельные значения обеих коэффициентов при к ^ 0 имеют 6 одинаковых значащих цифр (при использовании для экстраполяции значений при N0=16, 18, 20 отличие предельных коэффициентов около 0.1%). Численные значения коэффициента к\ ближе к предельным значениям, хотя для них отклонение от теоремы Гаусса немного выше. При к = 1/24 относительная погрешность при Vg = 1 около 4.4%, а при Vg = 0 больше (4.88%). Эти факты свидетельствуют о малой погрешности отношения сил при к < 1/24 и Ь —10, а также о полезности экстраполяции результатов.
Основные результаты были получены при определении силы взаимодействия по формулам (1.15), (1.16). Для контроля были выполнены расчеты коэффициента к при вычислении силы взаимодействия через потенциальные энергии (2.3). Этот способ оказался более трудоемким и более грубым. Он требовал расчета при двух близких значениях расстояния между шайбами. Кроме того, использование квадратной сетки не позволяло использовать произвольные значения числа узлов на радиусе шайбы N0. Расчеты (Ь = 10 и Ь = 9.8 с N0 < 24) показали, что вычисленное значение коэффициента к находится в пределах от 1.06 до 1.13, что не противоречит определенному выше предельному значению этого коэффициента 1.12575.
Отметим, что контроль решения по выполнимости теоремы Гаусса не дал полностью удовлетворительных результатов при равенстве зарядов (в случае + ф2 = 0 выполнимость теоремы следует из симметрии задачи). Величина ОО , показывающая относительное отклонение от выполнимости теоремы Гаусса, монотонно уменьшалась при уменьшении шага пространственной сетки. Это хорошо. Однако аппроксимация этой зависимости полиномом второй степени не давала нулевых предельных значений. При Vg = 0 предельное значение ОО и 1.26%, а при Vg = 1 ОО и 1.83%.
Выше подробно были обсуждены результаты для случая, когда расстояние между шайбами
значительно превосходило их размеры ( Ь = 10 ). Аналогичная серия расчетов была выполнена для расстояния между ближайшими торцами шайбами при Ь = 1. Соответствующие результаты представлены в табл. 2.
Таблица 2. Зависимость коэффициента
к (N0) для Ь = 1
N 0 16 20 24
к(0) 1.917132 1.90807 1.899554
к(1) 1.832624 1.825352 1.818128
дк 0.084508 0.082718 0.081426
продолжение табл. 2
N 0 26 28 30 Иш к
к(0) 1.895585 1.891813 1.888232 1.80768
к(1) 1.814773 1.811502 1.808373 1.73589
дк 0.080812 0.080311 0.079899 0.07173
Как видно, положительная разность значений Дк1= к(0) -к(1) монотонно убывает с ростом N0 . Однако убыль этой разницы происходит значительно медленнее, чем при Ь = 10 . Для относительной разницы получена оценка 8кх < 90 • к %. В последнем столбце таблицы, как и ранее, указаны предельные значения коэффициентов, полученные с помощью аппроксимации полиномом второй степени, построенного по трем наименьшим значениям шага сетки. Из анализа данных табл. 2, следует, что при Ь = 1 относительная погрешность определения отношения сил к значительна (8кг и 4.5% при N0 = 24 ). Погрешность менее одного процента характерна лишь для долей заряда и коэффициента к. В случае Vg ф 0 для уточнения предельных значений были дополнительно получены решения при N0 = 28 и N0 = 30 . Эти результаты позволили выяснить, что предельные значения по данным для N0 = 24,26,28 и N0 = 26,28, 30 отличаются менее чем на 0.01%.
Большую погрешность численного решения следует ожидать для еще меньших значений Ь < 1 (наименьшее значение Ь в расчетах равнялось 0.25). Значительная погрешность при малых расстояниях между шайбами ожидаема, так как при малых Ь в интервал между шайбами попадает малое число узлов. Упомянем об одной особенности итерационного процесса при малых расстояниях между шайбами. Уменьшение Ь уменьшает число узлов пространственной сетки на величину ДЬ • л • N0, что обычно ведет к сокращению числа итераций. Однако уменьшение Ь сокращает связи области II (см. рис. 1) со всей областью. Это обстоятельство может привести к увеличению числа итераций. В результате вычислительных экспериментов
выяснилось, что существует расстояние L » 2, при котором число итераций минимально, а любое изменение L увеличивает число итераций.
Как уже отмечалось, рассматриваемая геометрия имеет прямые углы. Известно, что в угловых точках напряженность поля резко возрастает. Для оценки эффекта угловых точек, была составлена специальная программа. В ней прямой угол заменялся на окружность радиуса R = h. Счет по этой программе (при L = 10, N 0 = 10) позволил выяснить, что интегральные характеристики - k, Vs1 = Vs2 = W отличаются от предельных менее чем на 7%, а распределения зарядов Qx1, QR1, Qx2 изменились в пределах 1%. Сглаживание границ уменьшило отношение напряженностей поля на границе k3 = max En /minEn до величины 2.942 (было
3.623). Выяснено, что вклад в силу от кругового сегмента существенен (около 35%). Любопытно отметить, что погрешность выполнимости теоремы Гаусса изменила знак и по модулю стала меньше в 2 раза. Расчеты по этой программе с меньшими значениями к не выполнялись, так как радиус сглаживания менялся вместе с шагом (Я0 = к) и серия расчетов с различными к давала решения для разных геометрий.
Перейдем к обсуждению главной зависимости к (Ь). В результате многочисленных вычислительных экспериментов (ВЭ) были выбраны 15 значений Ь , для которых были
определены предельные значения коэффициента к1 . Результаты этих ВЭ для определения зависимости к (Ь) представлены в табл.3.
Таблица 3. Зависимость к (Ь)
i 1 2 3 4 5 6 7 8
LH 10 9.0 8 7 6 5 4 3
k[ ] 1.12175 1.12408 1.12749 1.13251 1.1407 1.15426 1.1790 1.2298
i 9 10 11 12 13 14 15
L[i ] 2.5 2.0 1.5 1.0 0.75 0.50 0.25
k[ ] 1.2766 1.3533 1.4886 1.7359 1.955 2.28 2.791
Для описания аналитической зависимости к1 (Ь) была выбрана функция в виде ряда с тремя параметрами:
к = / (Ь) = 1+ь /(Ь+1)+Ь /(Ь+1)2 + Ь /(Ь+1)3. (3.6)
Предел этой функции при Ь ^ ж , как и требуется, стремится к 1. Значение (Ь +1) соответствует расстоянию между центрами шайб при толщинах шайб, равных 1. Значения коэффициентов, найденных методом наименьших квадратов, таковы: Ьх = 0.8008, Ь2 = 0.4836, Ъъ = 1.68549 . Основной вклад в сумму ряда (3.6) дает слагаемое с коэффициентом Ь2 . Сумма отклонений 15 значений к от функции (3.6) равна 0.06396, а соответствующая сумма модулей отклонений в 5 раз больше. Максимальное отклонение табличного значения к1 от функции (3.6) 0.0437 соответствует наибольшему расстоянию, использованному в расчетах Ь = 10 .
Из зависимости (3.6) следует, что сила взаимодействия в два раз меньше значения силы, вычисленной по формуле Кулона ( к1 = 2 ) при расстоянии между ближайшими торцами шайб Ь и 0.685 + 0.01.
Зависимость (3.6) представлена на рис. 4 сплошной линией. Точками (точнее - кругами)
на рисунке указаны соответствующие значения коэффициента к из табл. 3. Пунктирная линия соответствует значению к =1. Как видно, функция (3.6) дает удовлетворительное согласие с вычисленными значениями коэффициента к . Наибольшее относительное отклонение функции (3.6) от вычисленных значений к наблюдается при расстояниях Ь >7.
Предел отношения сил к при Ь ^ 0 конечен и равен к = к(0) и 1 + Ь1 + Ь2 + Ь3 и 3.97. Это значение к соответствует тому, что при соприкосновении шайб ( Ь = 0 ) заряды на шайбах в основном располагаются на дальних торцах. Предельное значение к = к (0) было подтверждено счетом по специальной программе, позволяющей вести счет при Ь = 0 .
В дальнейшем предполагается выполнить расчеты конечно-разностным методом для проводящих тел в форме шаров. В отличие от цилиндров шары имеют лишь один геометрический параметр и не имеют угловых точек, являющихся концентраторами напряжений. Такие расчеты позволят провести сравнение с результатами [1], полученными методом электростатических изображений.
Рис. 4. Зависимость отношения сил кх от расстояния между шайбами
5. Выводы
1. Найдена зависимость отклонения силы взаимодействия заряженных проводников от формулы Кулона для точечных зарядов.
2. Получена оценка расстояния между ближайшими торцами шайб Ь и 0.685 + 0.01, при котором вычисленная сила взаимодействия в два раза отличается от предсказания по закону Кулона для точечных зарядов.
3. Предложен способ уточнения расчетов путем экстраполяции результатов на сравнительно грубых сетках.
4. Показано, что наиболее точные значения силы взаимодействия получаются при вычислении пондермоторной силы, а вычисление этой силы через производные от потенциальной энергии требует больших затрат.
Список литературы
1. Саранин В. А. О взаимодействии двух электрически заряженных проводящих шаров // УФН. 1999. Т. 169. С. 453-458.
2. Саранин В. А., Майер В.В. Теоретические и экспериментальные исследования взаимодействия двух проводящих заряженных шаров // УФН. 2010. Т. 180. С. 1109-1117.
3. Саранин В. А. Электростатические осцилляторы // УФН, т. 182, 2012. с. 747-758.
4. Гладкий С. Л., Тарунин Е. Л., Ясницкий Л.Н. Применение метода фиктивных канонических областей в задачах электростатики // Вестник Пермского университета, Сер.: Физика. 2011. Вып. 3 (18). С.96-102.
5. Детлаф А.А., Яворский Б.М. Курс физики. М.: ВШ, 2001. 384 с.
6. Самарский А.А. Теория разностных схем. М.: Наука, 1977. 656 с.
7. Тарунин Е.Л. Вычислительный эксперимент в задачах свободной конвекции / Иркутский ун-т, Иркутск. 1990. 228 с.
8. Ильин В.П. Численные методы решения задач электростатики // М.: ФМ, 1989. 336 с.
9. Джексон Дж. Классическая электродинамика // М.: Мир, 1965. 699 с.
Interaction of charged conducting bodies on small distances
E. L. Tarunin
Perm State University, Bukirev St., 15, 614990, Perm
A force of interaction of two charged conducting bodies was calculated by means of finite-difference method. The bodies had a form of shot cylinders (pucks). The problem has axis symmetry. A potential of the electric field was found from solving of the corresponded Dirichxlet problem. It was investigated mainly the case of interaction of bodes with charges of the same sign. It was shown that the force of interaction for the case of small distances strongly differs from the corresponded Cou-lumb force. It was found the formula that describes the force of interaction as function on distance between bodies.
Keywords: Coulom's law, finite-difference method, Dirichxlet problem.