ПРИЛОЖЕНИЕ: Входное сопротивление нессиметричной антенны в среде с потерями




Теоретический анализ локального ближнепольного теплового эффекта при СВЧ сверления

 

E.Jerby, O.Aktushev, V.Dichtyar

(Journal of Applied Physics, 2004, vol.97)

 

Принцип СВЧ сверления основан на эффекте локального точечного нагрева (hot-spot), индуцированного коаксиальным ближнепольным нагревателем. СВЧ сверло локально расплавляет неметаллический материал и механически проникает во внутрь, образуя отверстие. Представленная статья демонстрирует теоретический анализ эффекта теплового нагрева, индуцированного оконечностью СВЧ сверла. Модель объединяет уравнения Максвелла и теплопроводности, включая свойства тепловой зависимости материала. При этом используется алгоритм конечных разностей во временной области (FDTD), примененный в двойной шкале времени. Моделирование продемонстрировано на силикате алюминия (муллите) и протестировано в упрощенных случаях. Полученные результаты показывают возрастание температуры от 1000 до 1300 К/с в области нагрева шириной 4 мм (приблизительно 0,1 длины волны СВЧ излучения). Порт входа отвечает за эффект ближнего поля, который моделируется эквивалентной схемой с изменяющимися во времени сосредоточенными параметрами. Помимо физического представления, теоретическое исследование обеспечивает вычислительные средства для проектирования и анализа СВЧ дрелей, мониторинга в реальном времени и адаптивного соответствия импеданса.

 

ВВЕДЕНИЕ

 

Неуправляемый нагрев (тепловой пробой) и место локального перегрева могут возникнуть случайно в различных процессах объемного СВЧ нагрева [1-6]. Материалы, восприимчивые к этим эффектам, характеризуются частично температурно-зависимыми свойствами, такими как зависимости увеличения потерь в диэлектрике или уменьшения теплопроводности во время повышения температуры. Поглощение СВЧ мощности приводит к концентрации температуры в определенной области, которая быстро превращается в область разогрева. В этом локальном самовозрастающем объединенном процессе (аналогично позитивной нелинейной обратной связи) рост температуры в локальной области приводит к точке плавления. Этот эффект неуправляемого нагрева может привести к серьезным повреждениям материала во время его локального СВЧ нагрева таким как спекание поверхности или сушка [7,8]. В этой работе, наоборот, процессы неуправляемого нагрева и локального перегрева считаются полезными для устройств локального нагрева и локальной плавки.

Метод СВЧ сверления использует эффект теплового пробоя для сверления твердых неметаллических материалов. Принципиальная схема СВЧ дрели показана на рис.1 (без источника СВЧ и состыковочных компонентов). СВЧ дрель концентрирует излучение энергии в ближнем поле в пределах материала и генерирует малую область локального нагрева, намного меньше длина волны СВЧ излучения. В процессе теплового пробоя материал становиться мягче или плавится изнутри в пределах локальной области нагрева. Затем центральный электрод СВЧ нагревателя внедряется в область локального нагрева для образования отверстия.

 

Рисунок 1 – Упрощенная схема электрода СВЧ дрели, состоящая из коаксиального волновода, заканчивающимся центральным электродом – «концентратором». Последний играет роль ближнепольной нессиметричной (или однополярной) антенны, (h) излучающей СВЧ энергию в область локального нагрева, приводящую к процессу сверления. Центральный электрод впоследствии механически внедряется в размягченный материал для образования отверстия.

 

Экспериментальные результаты, представленные в [9-13], демонстрируют эффективность СВЧ дрели в различных неметаллических материалах, включая бетон, силикатные минералы, кремний, стекла и керамику. Диаметр отверстий варьируется в диапазоне от 0,5 мм в стекле до 12 мм в бетоне. Эффект СВЧ нагрева достигается при концентрации СВЧ мощности обычно приблизительно до 102 Вт/мм2 в течении секунды. К примеру, СВЧ сверло может сделать отверстие или щель (канавку) в диэлектрическом покрытии без повреждений нижележащей металлической основы. Это продемонстрировано на теплобарьерной пленке из диоксида циркония (температура плавления >1500°C). В данной работе проведен теоретический анализ процесса теплового пробоя при помощи СВЧ сверла, с использованием конечно разностного метода анализа во временной области (FDTD - finite-difference time-domain).

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

 

ТЕОРЕТИЧЕСКАЯ МОДЕЛЬ

 

Анализ механизма СВЧ сверла объединяет электромагнитные (ЭМ) волны с вытекающими термическими эффектами [14,15]. ЭМ волны, излучаемые СВЧ нессиметричной антенной и мощность, поглощенная материалом, описываются уравнениями Максвелла для среды с потерями. Для структуры с цилиндрической симметрией продемонстрированной на рис.1 эти уравнения записываются в следующем виде:

где E r и E z – радиальные и продольные компоненты электрического поля соответственно; H φ – азимутальный компонент магнитного поля, σ d= ω 0 ε 0 ε" проводимость материала, связанная с диэлектрическими потерями; ε = ε 0(ε' - ") – комплексная диэлектрическая проницаемость материала.

Плотность поглощенной СВЧ мощности приблизительно определяется:

где <∙∙∙> обозначает локальное среднеквадратичное отклонение, усредненное по времени.

Уравнение теплопроводности служит для нахождения локального температурного профиля, возбужденного поглощенной СВЧ энергией. Уравнение теплопроводности для этого случая записано в следующем виде:

где ρ m, c m, k t - плотность тока, удельная теплоемкость и теплопроводность, соответственно; T – изменяющейся во времени температурный профиль.

Температурные зависимости от свойств среды (ρ m(T), c m(T), k t(T)) включены в данный анализ. Дополнительный нелинейный член в правой части уравнения (3) связан с температурной зависимостью k t(T) и его последующей пространственной вариации. В качестве примера температурные свойства температурно-зависимого диэлектрика (муллита) [17] показаны на рис.2. Они используются в 3тей части статьи для демонстрации эффекта теплового пробоя на муллите с использованием СВЧ дрели. Изменение свойств материала из-за локального нагрева приводит к дальнейшему локальному нагреву.

Рисунок 2 – Свойства температурных зависимостей муллита [17]: действительной и мнимой частей диэлектрической проницаемости (а); (b) удельной теплоемкости и теплопроводности.

 

В работе используется объединенный тепловой электромагнитный алгоритм, состоящий из двух численных решений, одно для электромагнитных уравнений (1) и (2), и другое для тепловых уравнений (3) [14,15]. Временная шкала тепловых эффектов намного шире, чем та же для ЭМ волн, обычно на 106 или больше (ЭМ волны распространяются по все системе за несколько наносекунд). Следовательно, приближение в двукратном масштабе времени применено для упрощения решения набора связанных уравнений (1)-(3). Распространение ЭМ волны рассчитывается в течении относительно короткого периода времени для каждого цикла вычисления уравнения теплопроводности. Свойства материала считаются постоянными для этих периодов времени вычисления. Результирующая поглощенная энергия обеспечивает начальные данные для уравнения теплопроводности (3). Локальная температура и параметры материала обновляются для каждого цикла вычисления уравнения теплопроводности. Температура может возрастать до точки плавления. Эффект фазового перехода (плавления или испарения) не включены в представленную модель.

Численное решение осуществляется по FDTD методу. Вначале задаются геометрические и начальные условия. ЭМ уравнение вычисляет сначала компонент магнитного поля (1с) а потом, использует результат для расчета компоненты электрического поля (1а) и (1b). Затем поглощение мощности ЭМ (2) рассчитывается для того, чтобы обеспечить генерирующий электрический источник для решения уравнения теплопроводности (3). Локальное поглощение мощности зависит от значения диэлектрической проницаемости, которое в свою очередь зависит от температуры. Значения диэлектрической проницаемости обновляются в каждом цикле расчета, в каждой ячейке дискретной FDTD решетки, принимая во внимание действующую локальную температуру. В дискретной форме FDTD уравнение (1b) к примеру выглядит следующей форме:

где Δ t и n – шаг по времени и индекс соответственно; i и j – радиальные и осевые индексы; и Δ r и Δ z их соответственные размеры шагов.

Локальная плотность поглощенной мощности P d (2) вычисляется как входное значение уравнения теплопроводности в центре каждой дискретной ячейки (например, температурный узел). Средние квадратичные значения электрического поля в этих точках (узлах) найдены путем интерполяции. Дискретное FDTD уравнение (3) записано в следующей форме:

где α = τ HeatEM – отношение между периодами расчета уравнения теплопроводности и электродинамического уравнения в каждом цикле (здесь α~106). Если не указано другого ε', ε", ρ m, с m и k t получают свое действительное значение в уравнении (4) и (5) в каждом i, j – узле.

Рабочая поверхность СВЧ сверла рассеивает тепло на (а) конвекцию через поверхностные границы, (b) излучение с поверхности вблизи концентратора и (с) электропроводность центрального электрода. Рассеивающийся удельный тепловой поток находится по формуле q cond­= h f(T – Ta), где h f – конвективный коэффициент теплопередачи и T a – температура окружающей среды. Излучение тепла находится по формуле q Rad= σε s(T4 – Ta4), где σ – константа Стефана Больцмана и ε s – поверхностная эмиссия (излучательная способность). Оба компонента – конвективная теплопередача и потери на излучение вычитаются из P dn(i,j) (5) в соответствующей граничной ячейки. Граничные условия для уравнеий электродинамики могут быть выбраны или для идеального металла или для поглощающего граничного условия.

Вычисляются также ЭМ волны, распространяющиеся вдоль коаксиальной линии и питающие СВЧ сверло (рис 1). Результаты объединяют, преимущественно падающую и отраженную волны (отражения вызваны импедансным дисбалансом нагрузки СВЧ сверла). Диаграмма направленности результирующей стоячей волны позволяет найти коэффициент отражения, изменяющийся во времени Г D(t). Изменяющийся во времени импеданс СВЧ сверла ZD(t) находятся по:

где Z c=(2π)-1 Z 0ln(b / a) характеристическое сопротивление коаксиального волновода питающего СВЧ сверло; Z 001/201/2 – импеданс свободного пространства; a и b –внутренний и внешний радиусы коаксиальной структуры.

Эквивалентное нагрузочное сопротивление представленной в модели линии передачи, показанной на рис 3, раскладывается на реальную и мнимую составляющие (R D и X D, соответствено) в следующем виде:

где R heat – эффективная поглощенная мощность на тепло в области локального нагрева; R Rad – отвечает за неэффективную часть мощности, излученной вне области сверления. Реактивные элементы L и C представляют собой индуктивность и емкость, соответственно, в ближнепольной области СВЧ сверла.

Рисунок 3 – Эквивалентная схема СВЧ сверла. Область взаимодействия сверла и область сверлящегося материала представлены в виде эквивалентной схемы с сосредоточенными параметрами с изменяющимися во времени элементами, зависящими от процесса развития теплового нагрева (пробоя).

 

Эквивалентное нагрузочное сопротивление (7) изменяется во времени при увеличении температуры за счет выделяющегося тепла во время процесса СВЧ сверления. Отраженная СВЧ мощность пропорциональна |ГD(t)|2, следовательно требуется механизм стабилизации импеданса для того, чтобы оптимизировать эффективность СВЧ сверла.

В качестве показателей для FDTD моделирования, его результаты сравниваются в упрощённых случаях (например: неизменяющиеся или постоянные свойства материала) с аналитической электродинамической моделью и коммерческой программой моделирования. Аналитическая модель однополюсной (несимметричной) антенны в однородной среде с потерями была взята из [21,22], учитывая синусоидальный ток, текущий вдоль антенны. Используя теорему Пойнтинга, были найдены тепловые и излучательные компоненты, а также накопленная энергия однополюсной антенны. Результирующие аппроксимации для эквивалентных компонентов импеданса (7) представлены в приложении. FDTD. Результаты моделирования также сравнивались с результатами, полученными с использованием коммерческой программы ANSOFT-HFSS (версия 9).

Используя утилиты расчета, представленные выше, третий раздел посвящен анализу теплового процесса СВЧ сверления. Промоделирован переход температурного профиля от равномерного распространения к области локального разогрева и определён соответствующий импеданс СВЧ нагрузки.

 

 

АНАЛИЗ

 

Механизм СВЧ сверления и процесс ближнепольного теплового пробоя (перегрева) промоделированы в случае материала муллита. Выбранный материал – типичный пример керамик и силикатов с похожими свойствами. При анализе используются зависимости, представленные на рис.2. Другие параметры были выбраны следующим образом: ρ m = 2500 кг/м3, ε s=1 и h f =5 Вт/м2К. Температурно-зависимые параметры муллита (поликристаллической керамики) могут значительно варьироваться. Следовательно, данный анализ следует рассматривать как наглядную демонстрацию физических процессов, а не как точный количественный расчет для муллита в целом.

Размеры коаксиального волновода (рис.1) были выбраны следующие: внутренний диаметр центрального проводника коаксиала – 2 мм; внешний диаметр коаксиального волновода – 10 мм; длина коаксиала – 80 мм. Размеры радиусов коаксиала соответствуют характеристическому сопротивлению Z c≈100 Ом, а длина соответствует 2/3 длины волны на частоте 2,45 ГГц; входная СВЧ мощность – 800 Вт и температура окружающей среды T a=300K.

Температурный рост в передней части СВЧ сверла представлен на рис.4(а) для 600 и 800 Вт эффективной входной мощности. Результаты моделирования показывают различные коэффициенты нагрева для двух различных входных мощностей. При этом эффект теплового перегрева происходит только при 800 Вт входной мощности. При этом температура достигает ~900К. Скорость дальнейшего роста температуры возрастает от ~2∙102 К/с до 1,2∙103 К/с в течении двух секунд. При меньшей СВЧ мощности (600 Вт) наблюдается меньший рост температуры, при которой не достигается критическая точка теплового пробоя при тех же условиях. (Похожая граничная величина наблюдалась в других экспериментах СВЧ сверления при различных условиях [23]).

На рис.4(b) показывает ширину области теплового перегрева, как полную ширину на полумаксимуме (FWHM) температурного профиля. Примерно 4,5 мм ширина соответствует ~0,1 длины волны. Быстрый рост температуры в области ближнего поля связан с дальнейшим уменьшением области перегрева, как показывает кривая для 800 Вт на рис.4(b).Эта тенденция концентрировать энергию в области локального перегрева характеризует феномен теплового пробоя.

Рисунок 4 – Зависимости области локального перегрева под центральным электродом СВЧ сверла для 800Вт и 600Вт эффективной входной СВЧ мощности.

 

Температурный профиль, развившийся в ближнем поле СВЧ сверла, показан на рис.5(a). В этом примере центральный электрод введен на 5мм внутрь образца муллита. Результат модулирования после 2,6 секунд при 800Вт излучения показал концентрацию высокотемпературной области напротив центрального электрода. Соответствующие профили радиальных и осевых компонент электрического поля показаны на рис.5(b) и 5(c) соответственно. Максимальная напряженность электрического поля в области локального перегрева достигает значения 3∙105 В/м, что оказывается в три раза больше, чем максимум радиального электрического поля на стыке с коаксиальным волноводом. Аксиальный компонент электрического поля оказывается значительным из-за области перегрева, расположенной напротив центрального электрода; что естественно для правильной операции СВЧ сверления.

 

Рисунок 5 – Результаты модулирования температуры (а) и компонент поля (b,c) пространственного распределения в муллите. Центральный СВЧ электрод, введенный на глубину h=5мм (с), выглядит как черное пятно. Температурное распределение (а), радиальный и осевой компоненты ЭП показаны для 800 Вт мощности СВЧ излучения за 2,6 секунды

 

Влияние материалов, зависящих от температуры, на эволюцию теплового перегрева демонстрируется четырьмя искусственными материалами, которые идентичны муллиту кроме одного постоянного параметра в каждом материале – значения комнатной температуры. На рис.6 показаны кривые изменения температуры в этих вырожденных материалах в сравнении с муллитом для всех различных параметров, которые зависят от температуры. Результаты показывают позитивный вклад ε" (Т) в реальном муллите в эволюцию теплового перегрева, тогда как реальная c m(T) играет препятствующую роль в этом отношении. Эффект от влияния ε' (Т) и k t(T) менее значителен. Отметим что уменьшение влияния k t(T) обусловлен второй (нелинейной) частью в правой части уравнения (3). Эти результаты совпадают с нашим интуитивным физическим пониманием участия этих параметров в уравнениях (3) и (4). Таким образом, наблюдается повышение роста температуры как и тангенса угла диэлектрических потерь и уменьшение удельной теплоемкости (c m) с ростом температуры.

Рисунок 6 – Зависимость свойств материала от температуры, представленная для искусственных муллитоподобных материалов в сравнении с реальным муллитом. В каждой кривой такие параметры как cm0 сохранены постоянными с их значениями для комнатной температуры, тогда как другие параметры взяты из рис.2. Значения cm0, ε '0, k t0 и ε "0 были взяты при комнатной температуре для соответствующих параметров.

Эффективная входная мощность – 800 Вт.

 

Диапазон изменения нагрузочного импеданса СВЧ сверла с центральным электродом, проникающим на изменяющуюся глубину h, показан на рис.7(a) и 7(b) для реальной и мнимой части, соответственно. Результаты FDTD здесь сравнены с результатами, полученными с помощью HFSS и аналитической моделью (см. приложение) для профиля стандартной температуры (300К). Получено удовлетворительное сходство. (Отметим, что аналитическая модель в зачительно меньшей степени учитывает внешний радиус коаксиального волновода по сравнению с длинной нессиметричной антенны. Таким образом, мнимая часть откланяется на h<6мм). Глубина проникновения центрального электрода, для которой XD=0 (h≈11мм), обеспечивает резонансную длину антенны. Нагрузочный импеданс СВЧ сверла является либо емкостью (XD<0), либо индуктивностью (XD>0) для меньшей или большей глубины проникновения, соответственно. Импеданс СВЧ сверла в резонансе (RD≈14Ом при h≈11мм) значительно меньше, чем характеристическое сопротивление коаксиальной линии (Zc≈100Ом). Таким образом требуется четкое соответствие значения импеданса. Диапазон изменения RD (5-35 Ом), обусловленный изменяющейся длиной проникающего центрального электрода, может нуждаться в адаптивном механизме уравновешивания импеданса для того чтобы добиться оптимальной эффективности на протяжении всего процесса СВЧ сверления (Предположение фиксированной входной мощности в этом анализе считается с учетом изменяющегося согласования импедансов).

Рисунок 7 – Зависимость реальной (а) и мнимой (b) части импеданса СВЧ сверла от глубины проникновения центрального электрода: FDTD моделирование (■), проведенное при комнатной температуре (300К); HFSS моделирование (×); аналитическая модель (▲).

Быстрый рост температуры из за теплового перегрева возникает локально вблизи открытого конца коаксиальной линии. Так как это изменяет свойства материала, то это также может изменить импеданс нагрузки СВЧ сверла (около генератора) во время процесса нагрева. Резистивные компоненты эквивалентного сопротивления СВЧ сверла, выраженные как Rheat и RRad в уравнении (7), показаны на рис.8 для 2 мм фиксированной длины проникновения в муллит (см. рис 2). Различия между RHeat и RRad раскрывает действительное отношение между эффективной поглощенной мощностью (на нагрев) и потерями на излучение во время процесса СВЧ сверления. Так как область локального разогрева увеличивается, компонент эффективной мощности возрастает; таким образом больше энергии поглощается и преобразовывается в тепло. Эффект изменения импеданса из-за теплового перегрева на фиксированной длине проникновения менее важен, чем тот же эффект, вызванный изменяющейся длиной проникновения, что показано на рис.7(a) и 7(b). Реальный процесс СВЧ сверления [9-13] объединяет оба случая изменения импеданса.

Рисунок 8 – Зависимость эквивалентного сопротивления среды от временя для различных параметров (рис.2). Различия между RHeat и RRad показывают компоненты мощностей на эффективную часть и потери соответственно во время СВЧ сверления на фиксированной глубине проникновения h=2мм (при эффективной входной мощности 800 Вт).

Изменение коэффициента отражения ГD в зависимости от глубины проникновения электрода показано на рис.9. Две кривые показывают критические случаи при T=300K и T=1300K при 800 Вт эффективного входной мощности. Реальный коэффициент отражения изменяется по мере неоднородного изменения температуры от верхней кривой к нижней кривой. Коэффициент отражения ГD уменьшается, когда увеличивается глубина проникновения или когда увеличивается температура. Следовательно, когда происходит процесс СВЧ сверления, он приводит к улучшению согласования импедансов с точностью до коаксиальной линии (Частично потому, что увеличивается рост потерь на излучение). Таким образом этот псевдо эффект самосогласования был также рассмотрен в наших предварительных экспериментах с СВЧ сверлом.

ОБСУЖДЕНИЕ (ВЫВОДЫ)

Объединенная тепло электродинамическая FDTD модель СВЧ сверла была применена в нашей работе для анализа феномена СВЧ дрели на муллите (выбран для этого анализа как типичный представитель керамики и других похожих материалов). Результаты показывают развитие области локального разогрева напротив передней части СВЧ сверла, как процесс теплового пробоя. Были определены эффекты воздействия температурно-зависящих свойств материала на процесс разогрева. Также были найдены зависимости импеданса нагрузки СВЧ сверла от различной длины проникновения и изменения температуры. Анализ создает условия для проникновения в сущность физических принципов процесса СВЧ сверления.

Рисунок 9 – Зависимость коэффициента отражения ГD от глубины проникновения центрального электрода в муллит при температуре 300К (▲) и 1300К (♦) при 800 Вт эффективной входной мощности. Коэффициент отражения уменьшается при повышении температуры и по мере погружения электрода внутрь материала. Во время нестабильного изменения температуры реальное значение ГD(t) находится между этими кривыми; он начинается с верхней кривой при комнатной температуре и спускается вниз к нижней кривой во время увеличения локальной температуры, как показано на примере FDTD модулирования (●).

На практике эта модель может быть использована для проектирования различных механизмов СВЧ сверления для различных материалов и для анализа мониторинга их работы.

Дополнительные исследования, связанные с этим анализом также необходимы, например (a) численный анализ процесса СВЧ сверления для других материалов (учитывая, что информация из рисунка 2 может быть не доступна полностью); (b) экспериментальные исследования СВЧ сверления в различных материалах и условиях в сравнении с теорией; (c) дальнейшее развитие более общих моделей, описывающих весь процесс СВЧ сверления, включая эффекты фазового перехода (плавка и испарение) одновременно с проникновением центрального электрода в материал.

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

ПРИЛОЖЕНИЕ: Входное сопротивление нессиметричной антенны в среде с потерями

Это приложение представляет приближённые формулы для нагрузочного сопротивления нессиметричной антенны, помещенной в среду с потерями. Модель учитывает синусоидальный ток антенны и однородную среду [21,22]. Используя теорему Поинтинга и метод индуцированной ЭДС, примененных для рассеивающей среды; резистивные и реактивные компоненты импеданса (7) приблизительно равны:

где α = k (((1+ x 2)1/2-1)∙ ε' /2 )1/2 и β = k ∙(((1+ x 2)1/2+1)∙ ε' /2 )1/2 – коэффициент затухания и волновое число в среде, соответственно; x =tan δ = ε" / ε' – тангенс угла потерь; λ и k=2π/λ – длина волны в вакууме и волновое число, соответственно, и ε= ε' - ε" = ε' (1-j∙tan δ) – относительная диэлектрическая проницаемость среды. Эквивалентная длина антенны - Leq≈βh{1+0.19/[ln(h/a)-81]}, и Ka=60[ln(h/a)-1] - «Среднее» значение характеристического импеданса в свободном пространстве. Другие параметры – A=tanh[(α/β)∙Leq] и B=tanLeq. Модель учитывает, что размеры антенны удовлетворяют следующим условиям: h/a »10, b/a<2 и L=βh«π/2.

 

Список используемой литературы:



Поделиться:




Поиск по сайту

©2015-2024 poisk-ru.ru
Все права принадлежать их авторам. Данный сайт не претендует на авторства, а предоставляет бесплатное использование.
Дата создания страницы: 2016-04-26 Нарушение авторских прав и Нарушение персональных данных


Поиск по сайту: