РАДИОЭЛЕКТРОНИКА
УДК 517.946:681.3
И. А. Конников
МОДЕЛИРОВАНИЕ ЭЛЕКТРОМАГНИТНЫХ НАВОДОК В САПР ЭЛЕКТРОННЫХ МОДУЛЕЙ
Изложен подход к использованию классических методов радиотехники для моделирования помехонесущего электромагнитного поля в микроэлектронике, ориентированный на использование в САПР. Метод позволяет избежать интегрирования в комплексной области и снизить объем вычислений по сравнению с известными методами. Предложена простая эквивалентная схема, позволяющая точно воспроизвести время задержки переходного процесса при моделировании канала паразитной связи.
В настоящее время учет паразитных электромагнитных эффектов (ПЭМЭ) при разработке современных микроэлектронных модулей является "узким местом" систем автоматизированного проектирования (САПР) [1], ограничивающим размерность решаемых на надлежащем научно-техническом уровне проектных задач. Работы [1,2] посвящены решению указанной проблемы. В работе [1] обоснована необходимость прямого использования методов теории электромагнитного поля, которые органично учитывают распределенный характер системы, состоящей из канала паразитной связи, источника и приемника излучения, и рассмотрена специфика использования этих методов в решаемой задаче.
Предлагаемый подход. При решении поставленной задачи будем использовать методы математической физики и электротехники, адаптированные к специфике САПР микроэлектронных модулей. Решение ориентировано на проектные задачи большой размерности с максимальным использованием аналитических методов, реализуемых заранее при разработке математического обеспечения САПР, в отличие от численных методов, предполагающие проведение основного, причем гораздо большего, объема вычислений в процессе проектирования.
Решение волнового уравнения относительно поляризационного потенциала П элементарного диполя в ^-м слое слоистой среды давно получено и, как известно из работы [3], в случае осевой симметрии задачи описывается интегралом Зоммерфельда:
сю
П(г, г, ¿0) = М У М\г)Фи(А, г, ¿о)с1А, (1)
о
где г — радиус в цилиндрической системе координат; М — амплитудный множитель; 70 — функция Бесселя первого рода нулевого порядка; Л — параметр разделения [5], именуемый в работе [3, с. 504] параметром разложения; функция Ф^(Л,г,г0) — математическая модель слоистой среды [4], определяемая из граничных условий для поляризационного потенциала на границах раздела слоев [3, с. 503].
Проблема вычисления интегралов Зоммерфельда имеет долгую историю. Уже несколько поколений исследователей [2,3,5,6] отмечают, что интегральное представление решения (1) очень сложно для его практического использования. Известен ряд методов [2,3,6] приближенного вычисления интеграла (1), требующих достаточно высокой квалификации исследователей, причем получаемые результаты нередко имеют весьма ограниченную область корректного использования [3, с. 513] и требуют непомерно больших затрат машинного времени. Описанные трудности при расчете наводимых помех можно значительно уменьшить, используя эквивалентную постоянную распространения кэ, которая рассчитывается в квазистационарном приближении через отношение потенциалов, создаваемых источником в слоистой среде и свободном пространстве, и позволяет учитывать влияние неоднородности (слоистости) среды. Такое приближение достаточно корректно, если расстояние между границами раздела крайних слоев, формирующих учитываемую при анализе среду, значительно меньше длины волны, так как при \кЯ\ ^ 1 (к — постоянная распространения среды; Я — расстояние между элементарным источником поля и точкой, где вычисляется поле) превалирует статическая составляющая поля, обратно пропорциональная кубу расстояния Я. Функции Грина для электростатической и магнитостатической задач вычисляются достаточно несложно по известной методике [4,7, 8], разработанной на основе решения уравнения Лапласа. При таком подходе нет необходимости в интегрировании по переменной Л в комплексной области и вычисления значительно упрощаются. В этом и заключается основная идея предлагаемого метода.
Для вычисления электродвижущей силы (ЭДС) помехи полученные аналитические выражения для помехонесущего поля подлежат интегрированию по объему приемника помехи, для чего целесообразно использовать хорошо верифицированнный метод тонкого провода, давно и широко применяемый в инженерной практике [9-11]. Тогда наведенные ЭДС вычисляют, используя характеристики невозмущенного поля, которое смоделировано при отсутствии объекта-приемника так же, как это делается при расчете паразитных емкостей и индуктив-ностей для моделирования ПЭМЭ. Такое допущение достаточно корректно, поскольку типичный диаметр цилиндра с вертикальным током
в микроэлектронике на два-пять десятичных порядков меньше1 длины волны, а расстояния между источником и приемником поля обычно намного больше их диаметров. Отказ от приближения тонкого провода приводит к необходимости решать задачи дифракции, что совершенно неприемлемо, поскольку сколько-нибудь строгое решение таких задач даже для наиболее простых слоистых сред приводит к результатам, которые пригодны лишь для проведения научных исследований и которые невозможно использовать в САПР микроэлектронных модулей в силу сложности этих результатов и большой размерности проектных задач [1].
Реализацию описанной идеи проиллюстрируем конкретным примером.
Считаем, что источник и приемник помехи представляют собой вертикальные (т.е. ориентированные параллельно оси аппликат и перпендикулярно границам раздела слоев) отрезки бесконечно тонкого провода. Аппликатами начала и конца провода-источника длиной 1и являются ги и ги + /и; гп и гп + 1п — то же для приемника помехи длиной 1п. Источник и приемник помехи полагаем электрически короткими, наведенные ЭДС вычисляются в приближении тонкого провода. Расстояние между источником и приемником в азимутальной плоскости обозначим га.
Пусть, как предложено ранее, локальная эквивалентная относительная диэлектрическая проницаемость
е е ^^
£у = 1 . .
2е 2е 0
где д£(А) — полученная при решении электростатической задачи математическая модель той же среды, для которой получена модель Фи.
Аналогично локальная эквивалентная относительная магнитная проницаемость будет иметь вид
£y = -R0j dz0 Jo(\r)q£(\)d\,
ß-э = Jo(\r)q^(\)d,\ I —,
где ^(А) — полученная при решении магнитостатической задачи математическая модель той же среды, для которой получена модель Фи.
Локальная постоянная распространения эквивалентной однородной среды кл = показывает, как изменяется набег фазы
1Так, например, луженый вывод микросхемы имеет диаметр приблизительно 0,6... 0,7 мм, перемычка межслойной коммутации в микросхеме имеет диаметр на два-три десятичных порядка меньше.
z
z
и
и
в рассматриваемой слоистой среде по сравнению со свободным пространством:
кл = ш
\
^и + 1и СЮ ^и+1и СЮ
£0^0 / dzo / Jo(Ar)q^(A)dA / / dzo / Jo(Ar)qe(A)dA.
Функция Грина, которая является решением волнового уравнения, для однородной среды со свойствами V-го слоя имеет вид
ехр(гА^ Л)
^в - —-• (2)
Для упрощения последующих расчетов распространение волны по прямым от каждого элементарного источника поля (точки с аппликатой го) к проводу-приемнику в точки с координатами {га,гп} и {га,гп + + /п} целесообразно описывать аппроксимирующим выражением вида (2) для однородной среды. Для этого необходимо характеризовать слоистую среду неизменным (интегральным) значением постоянной распространения кэ. Потребуем, чтобы постоянная распространения кэ обеспечивала на поверхности приемника такой же набег фазы, как и локальная постоянная кл, зависящая от координат. Разделим расстояние Л0п между элементарным источником поля (точкой с аппликатой ¿о) и точкой на поверхности приемника (точкой с аппликатой гп) на малые отрезки АЛ. Набег фазы в среде с постоянной кл на расстоянии АЛ равен клАЛ. Устремив максимальный из отрезков АЛ к нулю и перейдя к пределу, получим, что набег фазы на расстоянии Л0п равен линейному интегралу вдоль отрезка прямой, соединяющего названные точки:
R0n Га
O.R = кл(Л)^Л = ——
\
Zи+lи СЮ
£o^o / dzo J Jo(Ar)q^(A)dA -г-o-dr, (3)
2и+(и СЮ
/ dzo / Jo(Ar)qe(A)dA
2и o
где Л — у/г2 + (г - ¿о)2; Лоп - у7Г2 + (¿п - ¿о)2; г - го+фп-го)/га.
Отсюда постоянная кэ — «д/Л0п (от го она не зависит). Более подробно расчет эквивалентной постоянной распространения кэ рассмотрен ниже.
Значение потенциала <^п, создаваемого проводником-источником помехи, вычисляется интегрированием функции Грина по длине проводника с весом п(го), описывающим распределение зарядов по его длине. Это распределение в квазистационарном приближении равномерно [12]. Используя закон сохранения заряда в интегральной форме
o
o
z
z
и
и
а
o
o
[13, с. 377], получаем, что
п(го) = г//(ш/и),
где г = \/ — 1; I — ток в проводнике-источнике помехи; ш — круговая частота.
Поэтому на азимутальном расстоянии г в точке с аппликатой г
/г/ г ехр(гк ./Л) Св(^о)п(^о )^о = --Г -^— ¿¿о.
Заменив экспоненту первыми четырьмя членами ряда Маклоре-на (что вполне достаточно для электрически коротких проводников и малых расстояний микроэлектроники) и проинтегрировав, получим
^п(г,г) = г/[/о(/и,^и) + гкэЛ(/и ,^и) —
— кэ2 /2(/и,^и) — гкэ3/з(/и,^и)^ (4п^ ш/и), (4)
где функции /о(/и,^и) = АгеЪ¿и + ^—" + АгеЪ^ г^и; /1(/и,ги) = /и; /2(/и,^и) = [г2/о + (г — ¿и\/ г2 + (г — ¿и)2 +
+ (¿и + /и — гVг2 + (¿и + /и — ¿)2 ]/4; /з(/и,ги) = г2/и/6 + [(г — ¿и)3 + + (ги + /и — г)3]/18.
При г = га в приближении тонкого провода создающая ток помехи разность потенциалов Дип = ^п(га,гп + /п) — ^п(га,гп).
Задача решена в общем виде, однако описанная методика имеет значительный резерв повышения точности по сравнению с изложенным иллюстрационным вариантом.
Использование разностной математической модели среды. Точность предлагаемой методики вычисления можно значительно повысить, если воспользоваться тем обстоятельством, что поле в слое, где расположен приемник помехи, определяется в основном физическими характеристиками этого слоя; поле вблизи границы раздела определяется в основном физическими характеристиками слоев, прилегающих к границе раздела. Влияние остальных слоев незначительно, особенно на малых расстояниях г, характерных для микроэлектроники. Для каналов связи, содержащих провод, этот феномен усиливается вследствие внешнего скин-эффекта [14, с. 259], т.е. эффекта концентрации внешнего поля около поверхности провода. С учетом этого поле в слоистой среде можно рассчитать, представив выражение (1) с помощью простейшего тождественного преобразования в виде главной составляющей, легко вычисляемой и делающей основной вклад в фор-
мирование поля помехи, и составляющей, вычисляемой приближенно:
сю
П(т, г, го) = Пг(т, г, г0) + М ^ Л(Ат)Пр(А, г, го)д,\ (5)
о
где П^г — главная составляющая; ПР(А, г, г0) = Ф^(А, г, г0) — f (А, г, г0) — разностная математическая модель среды, позволяющая учесть влияние факторов, не учтенных главной составляющей; f (А, г, г0) определяется аналогично Ф^ (А, г, г0) для физической модели среды, которая соответствует главной составляющей.
Используя формулу, приведенную в работе [3, с. 509], решение волнового уравнения для поляризационного потенциала поля у плоской границы раздела двух однородных полупространств несложно привести к виду:
П(т, г, г0) = П^г (т, г, г0) + --4 (т, г, ¿0),
где
с
4(т,г,г0) = 2 ! Л(Ат)ехр[ — (г + г + 0)^/А2 — к2 ]х 0
х (-, к2-1-, — , 1 ^А:
+ А2 — к2-1 уА—Щ ) '
п _ ildzo
4пшеь
exp(ikv Rv) + exp(ikv Rv-i)
— главная составляю-
Rv Rv-1
щая, которая описывает поле диполя и его зеркального изображения в однородной среде со свойствами верхнего полупространства; к = \/цеш2 + 1цош — постоянная распространения электромагнитной волны в среде, где вычисляется поле; е, ц, о — абсолютная диэлектрическая проницаемость, абсолютная магнитная проницаемость и удельная активная проводимость среды соответственно; Rv, Rv-1 — расстояния от точки, где вычисляется поле, до диполя и его зеркального изображения соответственно; д,г0 — длина элементарного диполя; индекс V — 1 при величинах R, е, ц, о, к указывает на их принадлежность нижнему полупространству (под границей раздела), индекс V — верхнему полупространству; несобственный интеграл 4 через известные функции в явном виде не выражается.
Значение Пг (значение главной составляющей поляризационного потенциала) легко вычисляется и в основном влияет на результат вычисления П, особенно при \к'^-1| ^ |к2|, что обычно имеет место на практике. Все остальные составляющие потенциала П объединены с 4 в один интеграл — второе слагаемое выражения (5), значение
которого (поправку) приходится вычислять приближенно, с помощью эквивалентной постоянной распространения кэ.
Теперь необходимо конкретизировать некоторые аспекты вычисления эквивалентной (интегральной) постоянной распространения слоистой среды.
Расчет эквивалентной постоянной распространения. Постоянная распространения кэ рассчитывается на основе вычисления инте-
СЮ
грала С = J 7о(Лг)д(Л, г, го)^Л, в котором функция д = д^ для числи-
о
теля в формуле (3) и д = де для знаменателя. Чтобы вычислить этот интеграл его следует представить в виде
С = С1 + С2, (6)
N А^ с
где С1 = / ^о(Лг)((Л,г,го)^Л; С2 = / 7о(Лг)д(Л, г, ¿о)^Л;
^=1
Л^ = /г; — нули функции Бесселя 7о; Ло = 0; N — целое число полуволн функции Бесселя, при данном г укладывающихся на интервале [0,Ли].
Значение Ли выбирается так, чтобы при любых г и го, соответствующих интервалу их значений в решаемой задаче, выполнялось условие
((Лм ,г,го) = (с , (7)
где = Нш д(Л, г,го); значения N и Ли для обеих функций д одинаковы.
На практике условие (7) можно выполнить с любой требуемой точностью. Тогда величину С2 можно выразить через в(Лиг) — значение в-функции [4] аргумента Лиг:
С2 = дсв(Лм г). (8)
Значение в-функции при таком значении аргумента не обязательно равно нулю, однако, в подавляющем большинстве случаев интегралом С2 можно пренебречь по двум причинам. Во-первых, при последующем интегрировании функции Грина по координатам в подынтегральном выражении для С2 (см. (6)) неизбежно появляется сомножитель 1/Лп, где п — кратность интегрирования2, и при больших Л подынтегральная функция быстро убывает вместе с С2 по абсолютной величине. Во-вторых, при больших значениях аргумента Ли г значение
2Так, например, после интегрирования по г и в первообразной функции должен появиться сомножитель 1 /А2.
в-функции г) невелико [4] и для слоистых сред с плоскопараллельными и неограниченными в азимутальном направлении слоями обычно составляет3 = 0.
Впрочем полагать С2 = 0 совсем не обязательно. Поскольку соотношение (7) можно обеспечить с любой требуемой точностью, значение С2 легко вычисляется по формуле (8) через асимптотическое значение д^ функции д(Л) и значение в-функции [4].
Учитывая специфику функции Бесселя, для вычисления С1 в выражении (6) интервал интегрирования целесообразно разбить на N шагов и на каждом шаге использовать квадратурную формулу Лобат-то [15, с. 258] с тремя узлами, два из которых расположены на границах шага интегрирования. Алгебраическая степень точности такой формулы равна 3, что вполне достаточно. Поскольку границы шагов интегрирования совпадают с нулями функции Бесселя, отличным от нуля и подлежащим учету на каждом шаге будет из трех лишь одно слагаемое используемой квадратурной формулы, что значительно сокращает вычисления. Выражение для С1 в формуле (6) примет вид
л N
4 ^
= - Л^2^0(Ли2г)д(Ли2,г,г0),
3 V=1
где Лv2 = ^ - Лv-l)/2; Л12Г = 1,202412779 [16, с. 227].
При вычислении поправки А^п(гп,га) в точке с аппликатой г = гп формула (3) примет вид:
Га
шЯоп Г
®я =-
Га J 0
\
zn +1и N
Z V=1
2и -dr. (9)
ZH+1H N
/ Ay2 J0(Av2r)q£(Ay2;zn;z0)dz0
" V=1 ^и
Дробь под знаком радикала в формуле (9) особенностей не имеет и интегрирование по переменной г несложно выполнить при помощи одной из известных приближенных квадратурных формул, которые точны для алгебраических многочленов третьей степени, в том числе с помощью формулы Гаусса [15] или формулы Лобатто.
3 Строгое теоретическое доказательство этого соотношения для общего случая не проводилось, однако, методами вычислительного эксперимента удалось обнаружить лишь один класс функций г^, для которых = 0, а именно — класс функций для вычисления полей в плоскости источника, параллельной границам раздела слоев (х = х0). В этом случае поле вычисляется по методике [4] через нули ©-функции. При = 0 применение предлагаемой методики также возможно, однако ее экономичность несколько снижается.
Для малых расстояний г произведение Л^г невелико, условие (7) выполняется при N = 1 и формула (9) существенно упрощается:
= wR
0п 0п
q£ (Ai2,Zn,zo)dzo. (10)
При N = 1 как числитель, так и знаменатель будут иметь завышенные значения. Поскольку знаки погрешностей одинаковы, при вычислении дроби эти погрешности отчасти взаимно компенсируются, причем степень компенсации тем выше, чем ближе свойства среды к свойствам свободного пространства.
Аналогично при вычислении поправки + 1п,га) в точке с
аппликатой г = + /п на расстоянии Л0п = \/г2 + (гп + 1п — г0)2 по формуле (3)
=
=
0п
\
^0^0 J q^(A12)Zn+1njz0)dz01 J
qe(A12)Zn+1n)Z0)dz0-
(11)
Отсюда для поля, распространяющегося на расстояние Д0п, эквивалентная постоянная распространения кэ = ад/Д0п. Для лучей, сходящихся на другом конце проводника-приемника (или электрически короткого отрезка этого проводника), к = а'л/Л0п. Из выражений (10) и (11) нетрудно заметить, что обе постоянные кэ и к от аппликаты г0 не зависят. В случае цилиндрической или плоской волны кэ = аг/га.
При использовании метода разностной математической модели среды изложенная методика относится к вычислению поправки. Рассмотрим вычисление главной составляющей помехи.
Поляризационный потенциал поля связан со скалярным потенциалом известным соотношением [5, с. 444]. Главная составляющая скалярного потенциала поля в V-м слое вычисляется интегрированием по объему объекта-источника поля:
) = - П(хо ,yo)ds
divn
vr;
где $и — площадь поперечного сечения источника; п(х0,у0) — плотность распределения заряда по площади ; хо и уо — абсцисса и ордината источника.
z
z
и
и
S
и
В случае, когда диаметр4 объекта-источника мал, по сравнению с его длиной и расстоянием между источником и приемником га, практически без потери точности источник можно считать нитевидным и тогда
(г, г) = — J Пг (Б.
Zи
При формировании разностной математической модели среды было принято, что главная составляющая поляризационного потенциала Пг в формуле (5) представляет собой сумму двух экспонент. Поскольку Пг имеет только одну компоненту, ориентированную вдоль оси аппликат, то, выполнив дифференциальную операцию векторного анализа, заменив указанные экспоненты первыми четырьмя членами ряда Маклорена и проинтегрировав по г0, после приведения подобных членов для нитевидного источника получим:
И
уг(г,г) = --{и(1и,ги,г, г) — /4 Ни, —ги,г,г) — Ь(г,г) +
4пшеи
+ к2 [и(1и,ги,г,х) — 1е(—1и, —ги,г,г)] — 1к117(г,г ) —
— к^[18(1и,ги ,г,г) — Те Ни, —2и,г,г)]}. (12)
Задача несколько усложняется, если поперечный размер источника имеет тот же порядок, что и расстояние га (как это может иметь место в случае перемычек межслойной коммутации или при вычислении коэффициента передачи канала связи): тогда правая часть выражения (12) должна быть проинтегрирована по площади поперечного сечения Б, что, впрочем, не может вызвать принципиальных затруднений. Если же источник и приемник поля имеют одинаковые поперечные размеры £ (в направлении, перпендикулярном направлению распространения волны), то для учета поперечного размера канала связи £ и усреднения потенциала вдоль этого же направления двукратное интегрирование можно провести, используя метод средних геометрических расстояний (СГР), что позволит существенно упростить выражение для потенциала. При вычислении таких интегралов метод СГР обеспечивает высокую точность, если максимальный размер области, где вычисляется поле, превосходит интервал интегрирования хотя бы на десятичный порядок [17]; такое соотношение нередко выполняется. Примем такую ориентацию системы координат, при которой источник и приемник поля расположены на оси абсцисс. Тогда, проинтегрировав
Пг по г0 аналитически, а по ординатам у и у0 при помощи метода
4Объект-источник не обязательно является круглым и под диаметром понимается наибольшая из хорд его поперечного сечения.
СГР, получим, что потенциал поля помехи в V-м слое описывается формулой (12), в которой г = гх = -у/(ж — х0)2 + ¿2; tс = £ехр(—3/2) — СГР отрезка прямой, который имеет длину ¿, от самого себя.
Главная составляющая ЭДС помехи при сформулированных допущениях
иг = ^г(га,^п + /п) — ^Г (га,г„). (13)
Пример использования методики. Рассмотрим расчет ЭДС помехи, наведенной вертикальным проволочным выводом микросхемы, заключенной в металлический корпус, в таком же соседнем выводе. Соответствующая конструкция схематично приведена на рис. 1. Диэлектрики с диэлектрическими проницаемостями ^ и £2 полагаем непроводящими и немагнитными.
Вычисление главной составляющей ЭДС помехи иг проводится по формулам (12) и (13). Вычислим поправку. При расчете эквивалентной диэлектрической проницаемости £э в качестве главной составляющей примем потенциал поля диполя и его зеркального изображения в однородной среде со свойствами верхнего полупространства. Следуя рекомендациям работы [7], потенциал электрического поля в квазистационарном приближении для верхнего слоя структуры, представленной на рис. 1, будем описывать с помощью математической модели Ф^2 = ехр(—Л|г — г0|) + 82(Л)ехр[Л(г — г0)] + + р2(Л) ехр[Л(г0 — г)], а для нижнего слоя — с помощью модели Ф^1(Л,г,г0) = 81(Л)ехр[Л(г — г0)] + р1(Л) ехр[Л(г0 — г)]. Записав с помощью этих выражений граничные условия [7], получим систему интегральных уравнений относительно неизвестных (Л) и (Л) (V = 1, 2); число уравнений 4 необходимо и достаточно для отыскания всех неизвестных функций 8 и р. С помощью интеграла Фурье-Бесселя [5] переходим к системе линейных алгебраических уравнений;
Рис. 1. Физическая модель микросхемы на диэлектрической подложке в металлическом корпусе с вертикальными проволочными выводами
тогда граничные условия примут вид:
82(Л) ехр[Л(Н — 20)] + Р2(Л) ехр[Л(^0 — Н)] + ехр(Л|Н — 20|) = 0;
82(Л) ехр[Л(Л, — 20)] + Р2(Л) ехр[Л(^0 — + ехр(—Л|Л, — 20|) =
= 81 (Л) ехр[Л(Л, — 20)] + Р1(Л) ехр[Л(^0 — 81 (Л) ехр(—Л20) + Р1(Л) ехр(Л^0) = 0;
82(Л) ехр[Л(Л, — 20)] — Р2(Л) ехр[Л(^0 — + ехр(—Л|Л, — 20|)} =
= £!-(81 (Л) ехр[Л(Л, — 20)] — Р1(Л)ехр[Л(20 — Л.)]}. Решив эту систему уравнений, получим
2£2
82 =
(£1^ — £2){ехр[2Л(й — 20)] + ехр(—2Н)}'
1 — ехр[2Л(Н — й)] (14)
Р1 = —82 1 + ехр[2Л(Н — Ь)];
81 = —р1 ехр(—2Л20); Р2 = —82 ехр[2Л(Н — 20)] — 1,
„ 1 + ехр(—2^)1 — ехр[2Л(Н — й)1
где х/ =----—-
1 — ехр(—2^) 1+ ехр[2Л(Н — й)]'
Выделив с помощью тождества Вебера-Липшица [7] главные составляющие потенциалов, получим, что для расчета эквивалентной диэлектрической проницаемости £э надо использовать формальную разностную математическую модель среды де1(Л) = 81(Л) ехр[Л(2 — 20)] + + р1(Л) ехр[Л(20 — 2)] — ехр(—Л|2 + 20|) — ехр(—Л|2 — 20|) для нижнего слоя и модель де2(Л) = 82(Л)ехр[Л(2 — 20)] + р2(Л) ехр[Л(20 — 2)] — — ехр(—Л|2 + 2 + 0|) для верхнего слоя, которые вычисляются с помощью формул (14).
Формирование разностной математической модели конструкции для расчета эквивалентной магнитной проницаемости проводим аналогично. Поскольку оба слоя считаются немагнитными, при расчете используется только одна разностная математическая модель. Вектор-потенциал поля в квазистационарном приближении для обоих слоев структуры, представленной на рис.1, опишем функцией
Фа (Л, 2,20) = р(Л) ехр[Л(20 — 2)] + 8(Л) ехр[Л(2 — 20)] + ехр(—Л|2 — 201,
где
Р(Л) = ехр [—(Нр+—Л)Н— 1 — ехр(—ЛН);
s(A) = 1 - exp( -AH)
exp(-2Az0) — exp[A(H — z0)]
Выделив с помощью тождества Вебера-Липшица [7] главные составляющие вектор-потенциала, получим, что для расчета эквивалентной магнитной проницаемости надо использовать математическую модель
<?ДА) = Р(А) ехр[А(го - г)] + з(А) ехр[А(г - ¿о)] - ехр(-А|г + го\.
Поскольку диэлектрическая проницаемость меняется вдоль оси аппликат, с учетом обозначений координат на рис. 1 и формул (10) и (11) выражения для кэ и кЭ примут следующий вид:
кэ = ш
\
н
£o^o J q^(Ai2,H — l,Zo)d,Zo
H-l
h H '
j q£i(Ai2,H - l,Zo)dzo + ^ q£2(Ai2,H - l,Zo)dzo
H l h
(15)
kl = ш
\
н
£o V'0 q^(A12,H,Z0)dz0
H-l
H
(16)
qei(Ai2 , H, Zo)dZo + q£2(Ai2, H, Zo)dZo
Hl
Интегралы в формулах (15) и (16) легко вычислить с помощью квадратурной формулы Гаусса на два узла [15], после чего поправка А ип рассчитывается по формуле (4) как разность потенциалов концов
приемника:
Д-Un = ^n(Zn + 1п,Га) - Vn(Zn,Ta).
В результате ЭДС помехи ип = иг + Аип. В частотной области помеху можно моделировать, добавляя источник напряжения ип. Во временной области задача решается сложнее. Универсальным методом моделирования помехи как в частотной, так и во временной области будет использование эквивалентной электрической схемы канала связи.
Получение эквивалентной схемы. Параметры эквивалентной схемы беспроводного канала электромагнитной связи найдем, воспользовавшись известным методом Элмора [18] — методом моментов. Как показано в работе [1], этот метод описывает интегральные характеристики переходного процесса как в электрически коротких, так и в электрически длинных каналах связи. С формальной, математической точки зрения моменты импульсной характеристики аналогичны
h
h
Рис. 2. Эквивалентные схемы моделирующих четырехполюсников:
a — обобщенная; б — уточненная
коэффициентам разложения плотности распределения значений случайной величины в ряд Грама-Шарлье типа А [19], хорошо известный в математической статистике. Абсолютная величина первого момента по Элмору является широко принятой интегральной характеристикой переходного процесса — оценкой времени его задержки, что дает возможность потребовать от модели точного воспроизведения времени задержки. Моделировать канал связи будем с помощью двух П-образных четырехполюсников одинаковой морфологии (рис. 2, а). Первый из них должен моделировать распространение помехонесуще-го поля от точки с аппликатой 2и на поверхности источника (г = 0) в точку с аппликатой 2п на поверхности приемника (г = га), т.е. между крайними нижними точками канала связи. Второй четырехполюсник включается аналогично между крайними верхними точками источника и приемника. Для электрически коротких проводников этого должно быть достаточно, поскольку в квазистационарном приближении, на основе которого функционирует подсистема схемотехнического проектирования САПР, в соответствии с законом Ома для участка цепи ток помехи формируется разностью потенциалов его концов. При определении параметров моделирующих четырехполюсников ограничимся первыми двумя членами ряда Маклорена для коэффициента передачи канала связи, поскольку по Элмору члены с более высокими степенями частоты на время задержки не влияют [1, 18]. Принятое приближение тонкого провода позволяет считать, что наличие объекта-приемника не вызывает возмущения поля, поэтому можно полагать, что канал электромагнитной связи и моделирующие его четырехполюсники нагружены на характеристический импеданс среды который, как известно [3, с. 449], равен отношению напряженностей электрического (Еи) и магнитного (Н) полей. Аналитические выражения для Еи и Ни элементарного источника поля в однородной среде давно получены, хорошо известны [20, с. 244] и будут использованы при расчете импеданса
Рассчитаем импеданс 2х. В точке с аппликатой гп на поверхности приемника (г = Га) с учетом формул (4) и (12) потенциал
^Е(Га ,£п) = ^г(Га^п) + ^п (Га^п) = М [/4(^,^1) - ^(-¿и, -¿и) +
+ /о//и - гК/5 + гкэ], (17)
г/
где амплитудный множитель М =- [3, с. 514].
Создаваемая объектом-источником напряженность электрического поля в канале связи
E = —
grad+ / i(kvAr + k3An) , (18)
где АГ = гк^ПГ — вектор-потенциал элементарного источника поля, соответствующий главной составляющей поля; вычисление поправки Ап рассматривалось ранее (см. формулу (6) и далее).
По формулам (18) и (2) с учетом только членов, влияющих на формирование времени задержки, получим
Еи = [ - grad(- div Пг) - 1 , [ gradехр(гкэЯ)
zи zи
Напряженность магнитного поля Н связана с поляризационным потенциалом известным соотношением [5, с. 444]. Воспользовавшись для вычисления Еи и Н аналитическими выражениями для напряженно стей поля элементарного источника поля в однородной среде [20], выполнив дифференцирование по координатам, заменив экспоненты отрезками рядов Маклорена по степеням частоты, сохранив члены с показателями степени менее трех и проинтегрировав по х0 от до ги + 1и, а также введя обозначение р = гш, получим:
Zx = Еи/Н = (до + ра1)/(рЬ1 + р%), (19)
где ао = -0 + С2 - Сэ + (4; аг = £о£э^о^э/К;
Га
Ь1 = еУ\/ео£у(-(б + С7 + Се - Сэ); С1 = 2 + (—~-_ I )2 ;
Га + (^ ^и 1и)
Га Га Га
С . /"
2 ~~2 1 ; чэ —2 ; 7 ; | ТХо; ^4
rl + (z - Zи)2' Т2 + (z + Zи + ¡и)2' Т2 + (z + ZH)2'
Сб = \/r2 +(z - zu- ¡и)2 - Vr2 +(z - zu)2; Се = 2 + ( _и—¡ил;
Та + (z zu ¡и)
Z= z ZИ ; z z + ZИ + ¡и ; Z z + ZИ ;
Т2 + (z - Z^2' 8 Т2 + (z + Zи + lи)2' 9 Т2 + (z + Zи)2'
ZH+1H
62 = -w£vJ rot Gdz0; величина G определена формулой (6) при
q = q^.
Реализовав операцию rot и использовав приведенный вариант5 алгоритма вычисления G, получим
N ^и+^и
4 N Г
62 = -^vV^ A;;2Ji(Av2rW q^(Av2,z,zo) dzo.
3 V=1 ^
Функция q^ особенностей не имеет и интегрирование по переменной z0 несложно выполнить при помощи одной из известных приближенных квадратурных формул.
При получении формулы (19) учитывались только члены со степенями переменной p, участвующими в формировании времени задержки.
Для крайних нижних точек канала связи коэффициент передачи по напряжению
K (ш)
,z„) + г2 + ¿2^п)]/[^r(tc,Zn) + ^п (tcZn)],
где парциальные потенциалы и определяются по формулам (4) и (12), в которых в соответствии с методом СГР при вычислении потенциала приемника r = \/r2 + ¿2, а при вычислении потенциала источника r = tc.
Тогда с учетом формулы (12) имеет вид:
K (0) =
14(1и ,ZM,tc,ZM) I4( ^^ ZM,tc,ZM) + l0(tc,ZM)/lM
(20)
где
/о ,*.) = ArshZh + 2и - f + Arsh- ^ -
v/TFT^I v/TF+t!:
T {, \ Л 1
+/
Io(tc,ZH) = Arsh-—--h Arsh-
Время задержки переходного процесса каналом для крайних нижних точек
т =
14(1и)^и , /и, ^и^о^О + ^О^м^ОА
5 Выбор вариант алгоритма не принципиален и зависит от требуемой степени точности вычислений.
и
v^o[v^vVnUh{\/Т2 + ¿2, zn) + £э^э]
иЦи,^, /4( /и, -?к,\/Г2а ,гп) +1о(у/Га+|,г11)//и
(21)
Анализ соотношения (21) показывает, что время задержки по Эл-мору зависит как от свойств V-го слоя, так и от свойств слоистой среды в целом.
Для моделирования канала была принята симметричная П-образная схема. Как известно [21], коэффициент передачи по напряжению для такого четырехполюсника
Кп(р) = (1 + Уп%п + %п/2н )-1,
где сопротивление нагрузки Zн = Zx; остальные обозначения ясны из рис. 2 а.
Для получения параметров моделирующего четырехполюсника необходимо равенство времен задержки (по Элмору) (коэффициентов передачи на постоянном токе) для моделируемого канала и четырехполюсника. Как показывает анализ структуры выражений (20) и (21), требуемые равенства можно обеспечить простейшей схемой замещения, представленной на рис. 2, б. Тогда Уп = рС1, = 1/(рС2), где
С = Ь2 - а^/ао__.
С1 = таоК(0)/[1 - К(0)] - а1 Ь1/ао;
Ь2/а1 - Ь1/ао
С2 =
[К(0) - 1]/К(0) + тао/а1
Эти значения емкостей позволяют точно (по Элмору) смоделировать время задержки переходного процесса в канале электромагнитной связи, учитывая составлющие поля, в том числе поле излучения [20]. Параметры второго моделирующего четырехполюсника рассчитывают аналогично с той лишь разницей, что в расчетные формулы подставляют координаты крайних верхних точек канала связи. Использование предлагаемой схемы замещения канала связи позволяет существенно повысить размерность проектных задач микроэлектроники, решаемых на необходимом научно-техническом уровне.
Вывод. Описанный подход к вычислению интегралов Зоммер-фельда и предложенная эквивалентная схема могут использоваться для моделирования электромагнитного поля произвольно ориентированных излучателей в различных технических областях, в том числе при разработке математического и программного обеспечения САПР микроэлектронных модулей.
СПИСОК ЛИТЕРАТУРЫ
1. Конников И. А. Моделирование паразитных электромагнитных эффектов при автоматизированном проектировании электронных модулей // Информационные технологии. - 2007. - № 5. - С. 9-17.
2. Агапов С. В., Чермошенцев С. Ф. Методы и средства анализа и прогнозирования электромагнитных излучений от электронных средств // Информационные технологии. - 2003. - № 11. - С. 2-12.
3. Стрэттон Дж. А. Теория электромагнетизма. - М., Л.: ОГИЗ, 1948. -540 с.
4. К о н н и к о в И. А. Математическая модель конструкции микросхемы // Математическое моделирование. - 2007. - T. 19, № 4. - С. 37-44.
5. Тихонов А. Н., Самарский А. А. Уравнения математической физики. -М.: Наука, 1977.- 736 с.
6. К и н г Р., Смит Г. Антенны в материальных средах. - М.: Мир, 1984. -824 с.
7. К о н н и к о в И. А. Два способа вычисления функции Грина для уравнения Лапласа // Прикладная физика. - 2007. - № 2. - С. 17-24.
8. К о н н и к о в И. А. Индуктивность пленочных проводников в слоистых средах // - Судостроение. 1981. - № 11. - С. 27-28.
9. В э н с Э. Ф. Влияние электромагнитных помех на экранированные кабели. -М.: Радио и связь, 1982. - 120 с.
10. Д р а б к и н А. Л. Элементарные излучатели в проводящей среде // Радиотехника и электроника. - 1972. - № 2. - С. 264-267.
11. Лавров Г. А., Князев А. С. Приземные и подземные антенны. - М.: Сов. радио, 1965. - 260 с.
12. К о н н и к о в И. А. Влияние плотности распределения заряда на емкость прямоугольной пленки в слоистой среде // Электричество. - 2007. - № 3. -С. 37-41.
13. К а п л я н с к и й А. Е., Лысенко А. П., Полотовский Л. С. Теоретические основы электротехники. - М.: Высш. шк., 1972. - 448 с.
14. Зоммерфельд А. Электродинамика. - М.: ИЛ, 1958. - 502 с.
15. К р ы л о в В. И., Шульгина Л. Т. Справочная книга по численному интегрированию. - М.: Наука, 1966. - 372 с.
16. Справочник по специальным функциям / Под ред. М. Абрамовица, И. Стиган. - М.: Наука, 1979. - 832 с.
17. Конников И. А. Емкость тонкого проводника прямоугольного сечения // Авиакосмическое приборостроение. - 2006. - № 11. - С. 19-25.
18. E l m o r e W. C. The transient response of damped linear networks with particular regard to wideband amplifiers // Journal of Applied Physics. - 1948. - № 1. -P. 11-15.
19. Кендалл М., Стьюарт А. Теория распределений. - М.: Наука, 1966. - 588 с.
20. Пановский В., Филипс М. Классическая электродинамика. - М.: Физматгиз, 1963. - 432 с.
21. Б е с с о н о в Л. A. Линейные электрические цепи. 3-е изд., перераб. и доп. -М.: Высш. шк., 1983. - 336 с.
Статья поступила в редакцию 22.05.2007
Игорь Аркадьевич Конников родился в 1947 г., окончил Ленинградский институт авиационного приборостроения в 1972 г. Канд. техн. наук, старший преподаватель Санкт-Петербургского государственного университета культуры. Автор более 50 научных работ в области исследования электромагнитных эффектов в элементах микроэлектронной аппаратуры.
I.A. Konnikov (b. 1947) graduated from the Leningrad Institute for Aviation Instrument Engineering in 1972. Ph. D. (Eng.), senior teacher of the Saint-Petersburg State University for Culture. Author of more than 50 publications in the field of study of electromagnetic effects in elements of microelectronic apparatus.