Различия

Здесь показаны различия между выбранной ревизией и текущей версией данной страницы.

ru:tutorial [2012/07/22 22:04]
deinega [Кремниевая подложка]
ru:tutorial [2014/01/11 05:04] (текущий)
deinega [Наклонное падение]
Строка 1: Строка 1:
 +======Учебный курс======
 +
 +В этом разделе на тестовых примерах иллюстрируется, как с помощью EMTL проводить численные расчеты FDTD. Подробное описание используемых для этого классов находится в разделе [[documentation]]. 
 +/*Подробно о способах визуализации полученных результатов читайте в разделе [[visualisation]].*
 +======Единицы измерения======
 +
 +У нас в программе применяется запись уравнений Максвелла без коэффициентов типа <tex>\epsilon_0</tex>, <tex>\mu_0</tex>, <tex>c</tex>, <tex>4\pi</tex>. Это приводит к тому, что единицы измерения поля отличны от единиц измерения СИ, но никаких неудобств это не создает, поскольку все интересующие нас величины (такие, как коэффициент прохождения или сечение рассеяния) выражаются в виде отношения полей. 
 +
 +Скорость света полагается равной единице, а значит частота электромагнитной волны <tex>f=1/\lambda</tex> (например, волне с частотой 10 соответствует длина волны 0.1). Единицей измерения времени является время необходимое свету для того, чтобы пройти единицы длины.
 +
 +======Подключение EMTL======
 +
 +Подключение EMTL осуществляется с помощью заголовочного файла uiexp.h.
 +
 +<code cpp>
 +  #include "uiexp.h"
 +</code>  
 +
 +Перед использованием библиотеки вы должны провести ее внутреннюю инициализацию с помощью функции emInit, куда передаются аргументы командной строки.
 +
 +<code cpp>
 +  int main(int argc, char **argv){
 +    emInit(argc,argv);
 +    ...               // ваш код
 +    return 0;
 +  }
 +</code>  
 +
 +У функции emInit имеются третий и четвертый необязательные аргументы. В них можно указывать имена для директории, куда записываются текстовые выходные файлы, и для директории, куда записываются большие вспомогательные двоичные файлы. По умолчанию эти аргументы равны пустым строкам, а все файлы пишутся в ту директорию, где находится исполняемый файл.
 +
 +Проведение эксперимента FDTD осуществляется посредством класса uiExperiment, объявленного в файле uiexp.h. Нужно создать объект этого класса, а далее последовательно вызывать его функции, с помощью которых задается геометрия эксперимента, проводится расчет и последующий анализ.
 +
 +======Геометрические тела======
 +
 +Как создавать геометрические тела описано [[/ru/ivutils#geometricheskie_tela|здесь]].
 +
 +Если тела выходят за границы вычислительного объема, то рассматривается только та их часть, которая в него попадает.
 +
 +======Среды======
 +
 +Среде отвечает класс emMedium. В его конструктор передаются диэлектрическая проницаемость <tex> \epsilon </tex>, электропроводность <tex> \sigma </tex>, магнитная проницаемость <tex> \mu </tex> и магнитные потери <tex> \sigma *</tex>. По умолчанию их значения равны 1, 0, 1, 0. 
 +
 +FDTD не предполагает возможности табличного задания зависимости диэлектрической проницаемости от частоты. Однако в FDTD можно задать эту зависимость в виде произвольного числа членов в форме Друде и Лоренца:
 +<tex>
 +\varepsilon(\omega)=\varepsilon_{\infty}
 +-\sum_{p=1}^{N_D}\frac{\Delta \varepsilon_p \omega_p^2}{i\omega\gamma_p+\omega^2}
 ++\sum_{p=1}^{N_L}\frac{\Delta \varepsilon_p \omega_p^2}{\omega_p^2-2i\omega\gamma_p-\omega^2}
 +</tex>
 +
 +Например, золото можно задать следующим образом
 +
 +<code cpp>
 +  emMedium m(5.4);
 +  m.AddPolarizability(emDrude(.14*1e17,.103*1e+15/2,1));
 +  m.AddPolarizability(emLorentz(.427*1e+16,.87*1e+15,2.54*.268));
 +  m.AddPolarizability(emLorentz(.523*1e+16,1.32*1e+15,2.54*.732));
 +  m.calibrate(1/(299792458*1e+6));
 +</code>
 +
 +Параметры конструкторов для членов Друде и Лоренца такие: плазменная частота <tex> \omega_0</tex> и затухание <tex> \gamma</tex> и величина <tex>\Delta \varepsilon </tex>. Конкретные данные для золота взяты [[http://www.opticsinfobase.org/oe/abstract.cfm?uri=OE-15-26-18119|отсюда]].
 +
 +Задаваемая плазменная частота <tex> \omega_0</tex> и затухание <tex> \gamma</tex> должны измеряться в радианах на единицу измерения времени FDTD. В данном примере изначальные данные для золота имеют размерность радиан / секунду. Если единицей измерения длины FDTD является микрон, то,  соответственно, единицей измерения времени FDTD является время, за которое свет проходит один микрон. Для того, чтобы придать значениям <tex> \omega_0</tex>, <tex> \gamma</tex> требуемую размерность, мы воспользовались функцией calibrate, которая умножает <tex> \omega_0</tex> и <tex> \gamma</tex> (а также <tex>\sigma </tex> и <tex>\sigma *</tex>) на свой аргумент (в данном случае на обратную скорость света, сек/мкм).
 +
 +В случаях, когда у Вас имеется табличная зависимость диэлектрической проницаемости среды от частоты, ее необходимо аппроксимировать с помощью произвольного числа членов Друде, Лоренца. Для некоторых часто используемых сред (золото, серебро, кремний, вольфрам и т. д.) такая аппроксимация уже есть в файле opt_data.h. Для того, чтобы ей воспользоваться, нужно вызвать подходящую функцию, возвращающую объект emMedium, например
 +
 +<code cpp>
 +  emMedium au=GetAu(); // au - среда, отвечающая золоту
 +</code>
 +
 +вернет золото из предыдущего примера. В этом примере предполагается, что единица длины FDTD равна одному микрону. Если у вас единица длины FDTD равна 0.1 микрону, то нужно писать
 +
 +<code cpp>
 +  emMedium au=GetAu(0.1);
 +</code>
 +
 +В EMTL имеется также поддержка членов в форме Дебая, однако она сейчас недоступна.
 +
 +======Рассеяние на сфере======
 +
 +Создадим последовательно код, задающий эксперимент с рассеянием плоской волны на сфере. Задача о рассеянии на сфере имеет аналитическое решение, пользуясь которым, мы сможем удостовериться в том, что код работает правильно.
 +
 +Создадим объект этого класса uiExperiment:
 +
 +<code cpp>
 +  uiExperiment task;
 +</code>
 +
 +Вычислительным объемом у нас будет являться параллелепипед с противоположным координатами (-10,-10,-10) и (10,10,10).
 +
 +<code cpp>
 +  task.SetInternalSpace(Vector_3(-10),Vector_3(10));
 +</code>
 +
 +Зададим поглощающие граничные условия PML по всем трем направлениям.
 +
 +<code cpp>
 +  task.SetBC(BC_PML,BC_PML,BC_PML);
 +</code>
 +
 +Шаг по сетке dr будет равен единице.
 +
 +<code cpp>
 +  task.SetResolution(1);
 +</code>
 +
 +Второй необязательный аргумент функции SetResolution есть отношение шага по времени к шагу по сетке, значение по умолчанию для которого равно 0.5 (при этом значении FDTD работает устойчиво). Это означает, что шаг по времени в нашем случае равен <tex>dt=0.5dr=0.5 </tex>.
 +
 +Теперь мы можем сразу перейти к расчету:
 +
 +<code cpp>
 +  task.Calculate(200);
 +  task.Analyze();
 +</code>
 +
 +Аргумент функции Calculate равен времени FDTD, в течение которого будет выполняться эксперимент. Поскольку шаг по времени равен 0.5, количество временных итераций равно <tex>200/0.5=400</tex>.
 +
 +Итак, мы получили следующий код:
 +
 +<code cpp>
 +  #include "uiexp.h"
 +  int main(int argc, char **argv){
 +    emInit(argc,argv);
 +    uiExperiment task;
 +    task.SetInternalSpace(Vector_3(-10),Vector_3(10));
 +    task.SetBC(BC_PML,BC_PML,BC_PML);
 +    task.SetResolution(1);
 +    task.Calculate(200);
 +    return 0;
 +  }
 +</code>  
 +
 +После того, как мы скомпилируем и запустим его (о том, как это делается, написано в разделе "Установка"), программа выдаст сообщения, отражающие ход ее работы. Как будет из них видно, она честно выполнит свои 400 итераций (в сообщених выдается 399, потому что нумерация итераций начинается с нуля). 
 +
 +Как следует сообщений, размер сетки по каждому направлению равен 40 шагов. Это связано с тем, что к размеру внутреннего вычислительного объема, равного <tex>10 - (-10) = 20</tex>, прибавился еще задний и передний слои PML толщиной 10 шагов (10 шагов это оптимальное значение по умолчанию, которое можно менять с помощью функции ConfPML).
 +
 +Поскольку начальное значение поля на сетке равнялось нулю, а источников света в нашем эксперименте нет, то конечное значение поля тоже нуль. Фактически, мы считали вакуум, в котором ничего не происходит.
 +
 +Далее мы будем последовательно заполнять вычислительный объем телами, источниками и детекторами, вызывая различные функции класса uiExperiment. Создаваемый код нужно помещать перед вызовом функции Calculate.
 +
 +=====Total Field / Scattered Field=====
 +
 +Для генерации плоской волны мы будем используеть метод Total Field / Scattered Field (TF/SF). Область Total Field будет находится внутри параллелепипеда с координатами противоположных вершин (-7,-7,-7) и (7,7,7).
 +
 +<code cpp>
 +  task.AddTFSFBox(Vector_3(-7),Vector_3(7));
 +</code>
 +
 +В качестве падающей волны, используемой в области Total Field, будем использовать ограниченный по времени импульс с плоским фронтом, распространяющийся в направлении (0,0,1), и поляризованный вдоль направления (1,0,0):
 +
 +<code cpp>
 +  task.SetPlaneWave(Vector_3(0,0,1),Vector_3(1,0,0));
 +</code>
 +
 +Зададим временную зависимость падающего импульса в форме Беренгера:
 +
 +<code cpp>
 +  task.SetSignal(new Berenger(1,25,5,125));
 +</code>
 +
 +Форма Беренгера выглядит так:
 +
 +<tex>E(t) = E_0 (t-t_0) \exp \left[ -\left(\frac{t-t_0}{t_w} \right)^2 \right] </tex>, <tex>\left|t-t_0\right| \le T</tex>
 +
 +<tex>E(t) = 0</tex>, <tex> \left|t-t_0\right| > T </tex>
 +
 +В нашем примере <tex>E_0=1</tex>,<tex>t_0=25</tex>,<tex>t_w=5</tex>,<tex>T=125</tex>
 +
 +=====Детекторы=====
 +
 +Для того, чтобы измерять поля внутри вычислительного объема, нужны детекторы. Они считывают значения поля в тех точках, куда они помещены, никак при этом на него не влияя. Значения на детекторах получаются путем интерполяции по близлежащим узлам сетки. В ходе численного эксперимента детекторы записывают свои показания в двоичные файлы. Далее на основании этих двоичных файлов будут созданы читаемые текстовые файлы.
 +
 +Поместим наш первый детектор в центр вычислительного объема:
 +
 +<code cpp>
 +  task.AddDetectorSet("c",Vector_3(),Vector_3(),iVector_3(1));
 +</code>
 +
 +Первый параметр функции AddDetectorSet – имя множества детекторов, которое будет содержаться в названии соответствующих выходных файлов. Следующие три параметра задают трехмерную сетку, на которой расположены детекторы (левый нижний и правый верхний вершины сетки и количество сеточных шагов вдоль трех направлений). В нашем случае эта сетка вырождается в точку (0,0,0).
 +
 +Во время работы функции Calculate детектор будет записывать свои показания в двоичный файл binary.c.dat. Для того, чтобы из этого файла сделать текстовую таблицу, нужно после вызова Calculate вызвать функцию Analyze
 +
 +<code cpp>
 +  task.Analyze()
 +</code>
 +
 +В результате будет создан текстовый файл c_r0000.d, в котором будет записана история поля в детекторе в виде таблицы с колонками t, x, y, z, Ex, Ey, Ez, Hx, Hy, Hz для времени, координат и значений полей E и H.
 +
 +Для построения графиков на основании выходных текстовых файлов EMTL мы пользуемся популярной графической программой gnuplot, которая скачивается отсюда http://www.gnuplot.info.
 +
 +Например, для просмотра истории значения компоненты Ex во времени в файле c_r0000.d, нужно открыть гнуплот, с помощью команды cd перейти в нужную директорию (имя директории должно быть в кавычках), а потом в нем набрать
 +
 +  gnuplot> p 'c_r0000.d' u 1:5 w l
 +
 +В результате вы увидите заданный нами импульс Беренгера.
 +
 +{{tutorial_sphere:berenger.png}}
 +
 +У функции AddDetectorSet есть пятый необязательный параметр, куда можно помещать разные битовые флаги. 
 +
 +Например, для получения поля в частотном представлении нужно включить флаг DET_F:
 +
 +<code cpp>
 +   task.AddDetectorSet("cf",Vector_3(),Vector_3(),iVector_3(1),DET_F);
 +</code>
 +
 +В выходном текстовом файле cf_r0000.d будет представлена зависимость поля от частоты. Как видно из рисунка ниже, центральная частота нашего импульса примерно равна 0.04, что соответствует длине волны 25.
 +
 +{{tutorial_sphere:berenger_f.png}}
 +
 +=====Фильм с распространением волны=====
 +
 +Мы научились получать временную историю сигнала и частотную форму сигнала в точке. Для того, чтобы посмотреть поле во всем вычислительном объеме, можно покрыть его сеткой детекторов.
 +
 +Например, создадим двумерную сетку детекторов, параллельную плоскости yz и проходящую через центр координат (0,0,0).
 +
 +<code cpp>
 +  task.AddDetectorSet("f",Vector_3(-9.5),Vector_3(9.5),iVector_3(1,20,20),DET_TSTEP);
 +</code>
 +
 +Отметим, что поскольку количество узлов сетки вдоль ось x равно одному, x координата узлов равна (-9.5+9.5)/2=0. Если бы количество узлов было больше одного (как для случая других двум направлений), крайние узлы сетки оказались бы в точках -9.5 и 9.5. Если вам не очень понятно, что здесь имеется в виду, имеет смысл построить разные сетки детекторов и посмотреть на них в гнуплоте (о том, как это делать, будет идти речь ниже).
 +
 +По умолчанию на выходе программа создает текстовые файлы в каждом из которых находится история сигнала для определенного детектора. Для того, чтобы создать выходные файлы для каждого временного шага со значениями полей на всей сетке детекторов, нужно использовать флаг DET_TSTEP.
 +
 +В нашем случае это приведет к созданию файлов f_txxxx.d, где xxxx – номер временной итерации. Просмотрим распределение поля на 50 временной итерации, набрав в гнуплоте:
 +
 +  gnuplot> sp 'f_t0050.d' u 3:4:5
 +
 +{{tutorial_sphere:berenger_3d.png}}
 +
 +Для последовательного просмотра всех итераций создадим файл film.gnu и скопируем в него следующий скрипт
 +
 +  pref="f"
 +  fn=sprintf("\%s\_t\%04d.d",pref,i)
 +  splot [][][-3:3] fn u 3:4:5 w l
 +  i=i+1
 +  if(i<iterations)reread
 +
 +Запустим его, набрав в командной строке
 +
 +  gnuplot> i=0; iterations=400; load 'film.gnu'
 +
 +В итоге мы увидим плоскую волну в форме импульса Беренгера, проходящую в области Total Field.
 +
 +Фильм можно сделать более симпатичным, если значение поля Ex будет передаваться не высотой по оси z, а непрерывно меняющимся цветом. Для этого наберем в командной строке gnuplot
 +
 +  gnuplot> set pm3d interpolate 2,2
 +  gnuplot> cbr=3
 +  gnuplot> set palette defined (-cbr 'blue', 0 'white', cbr 'red')
 +  gnuplot> set cbrange [-cbr:cbr]
 +  gnuplot> set view 0,0
 +
 +а в скрипте film.gnu вместо
 +
 +  splot [][][-3:3] fn u 3:4:5 w l
 +
 +напишем
 +
 +  splot [][][-cbr:cbr] fn u 3:4:5 w pm3d
 +
 +Допустим, мы хотим создать видео. Для этого наберем в командной строке
 +
 +  gnuplot> set terminal png
 +
 +Теперь выходные данные gnuplotа будут выводиться не на экран, а в файл в формате 'png'.
 +
 +Модифицируем наш скрипт следующим образом
 +
 +  pref="f" 
 +  fp=sprintf("\%s\%04d.png",pref,i)
 +  set output fp
 +  fn=sprintf("\%s\_t\%04d.d",pref,i)
 +  splot [][][-cbr: cbr] fn u 3:4:5 w pm3d
 +  i=i+1
 +  if(i<iterations)reread
 +
 +Запустив его, мы получим 400 картинок, которые впоследствии можно слить в один файл avi (например, программой VirtualDub). 
 +
 +Чтобы вернуться в режим вывода на экран и вернуться из режима pm3d, нужно набрать
 +
 +  gnuplot> set terminal win
 +  gnuplot> unset pm3d
 +
 +Для того, чтобы получить подобные картинки в частотном представлении, нужно использовать сразу два битовых флага DET_F и DET_TSTEP:
 +
 +<code cpp>
 +  task.AddDetectorSet("ff",Vector_3(-9.5),Vector_3(9.5),iVector_3(1,20,20),DET_F|DET_TSTEP);
 +</code>
 +
 +Как видно, в С++ битовые флаги можно объединять с помощью оператора |.
 +
 +=====Добавление тела=====
 +
 +Добавим стеклянную (показатель преломления <tex>n=1.45</tex>, а диэлектрическая проницаемость <tex>\epsilon=n^2</tex>) сферу радиуса 5 в центр контейнера:
 +
 +<code cpp>
 +  task.AddObject(emMedium(1.45*1.45),new Sphere(5,Vector_3()));
 +</code>
 +
 +Желательно, чтобы между телом и границей TF/SF было расстояние не меньше 1.5 сеточных шагов по пространству. Мы это условие соблюдаем (7-5=2).
 +
 +После этого можно вновь создать фильм и увидеть, как сфера частично рассеивает падающую волну в область Scattered Field. Рассеянная волна через какое-то время поглощается PML, что свидетельствует о том, что эксперимент можно заканчивать. Что мы и делаем, выбрав время эксперимента 200.
 +
 +{{tutorial_sphere:film_little_mie.png?500}}
 +
 +Если бы вместо PML мы бы использовали отражающие граничные условия, то рассеянная волна переотражалась бы от границ контейнера.
 +
 +<code cpp>
 +  task.SetBC(BC_REFL,BC_REFL,BC_REFL);
 +</code>
 +
 +Это также можно увидеть с помощью фильма, который также можно построить в качестве домашнего задания.
 +
 +gnuplot позволяет рисовать линии, вставлять надписи, регулировать размеры осей, а также их наличие. Это позволяет делать достаточно красивые фильмы. Вот один из них.
 +
 +{{tutorial_sphere:mie_scattering.png?500}}
 +
 +В этом фильме свет падает под углом <tex>45^0</tex> к осям y,z. Для того, чтобы отчетливей была видна рассеянная сферой волна, мы увеличили расстояние между сферой, границей TF/SF и PML.
 +
 +=====Поток энергии через поверхность=====
 +
 +Прежде чем читать дальше, [[fdtd#Рассеяние на теле|ознакомьтесь с теорией]].
 +
 +Программа FDTD позволяет рассчитать поток энергии заданную поверхность. Для этого на этой поверхности размещается сетка детекторов, данные которых участвуют в расчете этого потока. Такой поверхностью может являться поверхность параллелепипеда
 +
 +<code cpp>
 +  task.AddFluxSet("flux",Vector_3(-8.5),Vector_3(8.5));
 +</code>
 +
 +Созданные детекторы будут находиться друг от друга примерно на расстоянии одного сеточного шага.
 +
 +Мы создали поверхность в области Scattered-Field для того, чтобы рассчитать поток энергии отраженной от сферы волны <tex>W_{sca}</tex>. Желательно, чтобы между детекторами и границей TF/SF было расстояние не меньше 1.5 сеточных шагов по пространству. Мы это условие соблюдаем (8.5-7=1.5).
 +
 +После анализа будет создан файл flux.d с двумя колонками для частоты и соответствующего потока энергии. Его можно использовать для расчета эффективности рассеяния <tex>Q_{sca}</tex> для нашей стеклянной сферы, который равен
 +
 +<tex>
 +Q_{sca}=\frac{W_{sca}}{I \dot G}
 +</tex>,
 +
 +где
 +  *<tex>W_{sca}</tex> - поток энергии рассеянной волны через окружающую сферу поверхность
 +  *<tex>I</tex> - интенсивность падающей волны
 +  *<tex>G = \pi 5^2</tex> - площадь проекции сферы на плоскость перпендикулярную направлению падающей волны.
 +
 +Для расчета интенсивности падающего света мы можем воспользоваться полученным ранее файлом cf_r0000.d, в котором хранится падающее поле в частотном представлении. Поскольку падающая плоская волна поляризована в направлении x, ее интенсивность можно получить как <tex>I=|\frac{1}{2}{E_x}^2|</tex>.
 +
 +Ниже придставлено сравнение для <tex>Q_{sca}</tex> между FDTD и программой BHMIE. Значения <tex>Q_{sca}</tex> представлены в логарифмическом масштабе.
 +
 +{{tutorial_sphere:mie_qsca.png?400}}
 +
 +Программа BHMIE основана на суммировании рядов, составляющих аналитическое решение задачи рассеяния на сфере (решение этой задачи было опубликовано Ми в 1908 г.). Описание программы BHMIE можно найти в книге ((C. F. Bohren and D. R. Huffman: Absorption and Scattering of Light by Small Particles, Wiley-Interscience, New York (1983) )). Сама программа (вместе с другими замечательными программами) скачивается отсюда: http://atol.ucsd.edu/scatlib/index.htm, где она любезно выложена Dr Piotr J. Flatau из University of California San Diego.
 +
 +=====Просмотр заданной геометрии=====
 +
 +Для того чтобы можно просматривать в gnuplot расположение задаваемых объектов, программа на выходе создает некоторые вспомогательные файлы. Перечислим их:
 +
 +  * cell.pol – вычислительный объем
 +  * med.pol – тела
 +  * pmli.pol, pmlo.pol - внутрення и наружняя границы PML 
 +  * tf.pol – область Total Field
 +  * sg.vec – направление падающей волны
 +  * файлы с суффиком _vset - детекторы
 +
 +Расширение pol в имени файлов означает, что в файлах задается многогранник(и) (polyhedron), расширение d  - точки (dot), а расширение vec - вектора (vector). 
 +
 +Наберем в командной строке gnuplotа
 +
 +  gnuplot> sp 'cell.pol' w l, 'med.pol' w l, 'pmlo.pol' w l, 'pmli.pol' w l, 'tf.pol' w l, 'sg.vec' w vec, 'c_vset.d', 'flux_vset.d' w d
 +
 +В результате мы увидим следующее.
 +
 +{{tutorial_sphere:container_mie3.png}}
 +
 +=====Преобразование ближнего поля в дальнее=====
 +
 +Прежде чем читать дальше, [[fdtd#Матрица амплитуд рассеяния|ознакомьтесь с теорией]].
 +
 +Для расчета углового распределения рассеянной волны используется матрица амплитуд рассеяния, связывающая амплитуду рассеянного поля в дальней зоне с амплитудой падающей волны. Для получения рассеянного поля в дальней зоне можно было бы значительно увеличить вычислительный объем, так чтобы он захватывал собой дальнюю зону. Однако имеется более экономный путь.
 +
 +Окружим тело мыслимой замкнутой поверхностью в области Scattered-Field. Значение поля вне этой поверхности в начальный момент времени <tex>t=0</tex> равно нулю. В течение численного эксперимента нам становится известно значение поля на поверхности при <tex>t\ge{}0</tex>. Из свойств уравнений Максвелла вытекает, что этого достаточно для однозначного определения поля вне поверхности при <tex>t\ge{}0</tex>. Это можно сделать с помощью техники преобразования ближнего поля в дальнее (Near-to-Far-Field Transformation). В EMTL эта техника реализована так, что на выходе вы получите дальние поля сразу в частотном представлении.
 +
 +Окружим тело замкнутой поверхностью в области Scattered-Field (она совпадает с поверхностью, использованной нами ранее для расчета потока рассеянной волны).
 +
 +<code cpp>
 +  Vector_3 origin; // центр, откуда отмеряется расстояние до "дальних" точек
 +  valtype length=1000; // расстояние до "дальних" точек
 +  // сетка по углам
 +  valtype th0=0,th1=M_PI;
 +  valtype fi0=0,fi1=M_PI/2;
 +  int thn=19,fin=2;
 +  task.AddNearToFarSet("far",Vector_3(-8.5),Vector_3(8.5),origin,length,th0,th1,thn,fi0,fi1,fin,DET_TSTEP);
 +</code>
 +
 +Значения полей на этой поверхности, полученные в ходе численного эксперимента, будут использоваться на этапе анализа для получения полей в точках, удаленных на расстояние 1000 от начала координат в направлениях, задаваемых углами <tex>\theta</tex> и <tex>\phi</tex>. Угол <tex>\theta</tex> пробегает значения от 0 до <tex>\pi</tex> с шагом <tex>\pi/18</tex>, а угол <tex>\phi</tex> принимает значения 0 и <tex>\pi/2</tex>.
 +
 +По окончанию анализа создадутся файлы 'far_xxxx.d, где xxxx – номер частоты (частоты нумерются по возрастанию). Поля <tex>\bf E</tex> и <tex>\bf H</tex> в этих файлах будут разложены по ортонормированному базису, связанному со сферической системой координат (см. [[fdtd#Матрица амплитуд рассеяния|описание]]).
 +
 +Сделаем несколько замечаний по поводу содержимого 'far.d'.
 +
 +Во-первых, <tex>E_r</tex> близко к нулю, что согласуется с тем фактом, что рассеянное поле в дальней зоне поперечно.
 +
 +Во-вторых, <tex>E_{\phi}</tex> близко к нулю для угла <tex>\phi=0</tex>, а <tex>E_{\theta}</tex> близко к нулю для угла <tex>\phi=\pi/2</tex>. Это связано с тем, что недиагональные элементы матрицы амплитуд рассеяния для сферы равны нулю: <tex>S_3=S_4=0</tex>.
 +
 +Отметим также, что при рассеянии на сфере матричные элементы <tex>S_1, S_2</tex> не зависят от угла <tex>\phi</tex>.
 +
 +Теперь, пользуясь абсолютным значением <tex>E_x</tex> из файла cf_r0000.d, а также принимая во внимание изложенную [[fdtd#Матрица амплитуд рассеяния|здесь]] теорию, мы можем вычислить абсолютные значения элементов <tex>\left|S_1\right|, \left|S_2\right|</tex>.
 +
 +Ниже представлено сравнение полученных результатов с результатами программы BHMIE для частоты <tex>f=0.05</tex>. Заметим, что для углов <tex>\theta</tex> близких к <tex>\pi</tex> согласие хуже. Это явление известно в литературе по FDTD под названием Backscattering Problem.
 +
 +{{tutorial_sphere:mie_c12.png?400}}
 +
 +=====Фазы работы=====
 +
 +Как можно видеть из вышесказанного, объект uiExperiment может последовательно находиться в следующих трех фазах:
 +
 +  * Построение геометрии. На этой фазе вы должны определить размер вычислительного объема, шаг по сетке и граничные условия, а также задать тела, источники и детекторы.
 +  * Расчет, протекающий во время работы функции Calculate
 +  * Анализ, вызываемый функцией Analyze. На этой фазе производится анализ выходных двоичных файлов, созданных в ходе расчета. На основании этих двоичных файлов создаются читаемые текстовые файлы.
 +
 +Иногда бывает полезно, чтобы были пройдены не все фазы. Для этого нужно использовать функцию SetPhases, параметр которой может принимать значения "g", "c" и "a". Например, если вы укажете
 +
 +<code cpp>
 +  task.SetPhases("g");
 +</code>
 +
 +то будет пройдена только фаза построения геометрии. Это полезно, если вы хотите ознакомиться с построенной геометрией в gnuplot без проведения самого расчета (например, чтобы проверить, правильно ли она задана).
 +
 +Значение "c" и приведут к тому, что будет выполнено построение геометрии и расчет. Это полезно, если вы хотите отложить анализ больших двоичных файлов на потом (например, когда задача будет перезапущена на другом числе процессоров).
 +
 +Значение "a" означает выполнение только фазы анализа.
 +
 +======Спектры прохождения и отражения для периодических структур======
 +
 +С помощью FDTD можно рассчитывать спектры прохождения и отражения для периодических структур (фотонные кристаллы, массивы антенн и т.д.)
 +Для их расчета рассчитывается распространение ограниченного по времени импульса через структуру.
 +Периодичность моделируется путем задания периодических граничных условий в нужных направлениях.
 +
 +Рассмотрим элементарную ячейку фотонно-кристаллической пластинки, состоящей из сфер упакованных в квадратную решетку.
 +Такая ячейка ограничена периодическими границами 1 и 2.
 +Граница Total Field / Scattered Field (TF/SF) 3 генерирует падаюший импульс (область Total Field заштрихована).
 +Для непериодического направления используются поглощающие граничные уловия PML.
 +Они поглощают отраженную и прошедшую волны, моделируя их уход на бесконечность.
 +Численный эксперимент ведется до тех пор пока излучение не выйдет из структуры.
 +Прошедшая и отраженная волны записываются во время эксперимента, преобразуются в частотное представление и нормализуются на спетр падающей волны, что дает спектры прохождения и отражения.
 +
 +{{tutorial_periodic:geometry_norm.png?250}}\\
 +
 +В случае нормального падения периодические граничные условия выглядят так
 +
 +<tex>
 +{\bf F}({\bf x},t)={\bf F}\left({\bf x} \pm {\bf a},t\right),
 +</tex>
 +
 +где <tex>\bf F</tex> это электрическое или магнитное поле (<tex>\bf E</tex> или <tex>\bf H</tex>), <tex>\bf a</tex> это вектор трансляции для периодической структуры, <tex>\bf x</tex> и <tex>t</tex> пространственная и временная координаты.
 +
 +В случае наклонного падения периодические граничные условия содержат временной сдвиг. 
 +В 2D они принимают форму
 +
 +<tex>
 +{\bf F}({\bf x},t)={\bf F}\left({\bf x}\pm{\bf a},t\pm a\sin{\theta}/c\right),
 +</tex>
 +
 +где <tex>\theta</tex> это угол падения, <tex>c</tex> - скорость света в среде, откуда падает свет.
 +Смысл этого выражения легко понять из следующих соображений.
 +Падающий импульс приходит на границу 1 раньше чем на границу 2, поэтому значения полей на этих границах "сдвинуты во времени" относительно друг друга.
 +Использование периодических граничных условий требует знания поля в прошлом на границе 1 (для применения на границе 2) и поля в будущем на границе 2 (для применения на границе 1).
 +Поле в прошлом могут быть получено, если ведется запись распространения волны.
 +Получение поля в будущем проблематично, поскольку оно неизвестно заранее.
 +Для решения этой проблемы можно использовать итерационный метод для моделирования наклонного падения на периодическую структуру.
 +
 +Ниже мы рассматрим моделирование нормального падения на одномерном примере, далее мы объясним принцир итерационного метода и покажем как он наботает на примере фотонно-кристаллической пластинки из сфер и антиотражающего текстурированного покрытия.
 +
 +
 +=====Бесконечная пластинка=====
 +
 +Ниже мы рассмотрим как рассчитать спектры прохождения и отражения для бесконечной пластинки (которая может рассматриваться как периодическая с периодом нуль) в случае нормального падения.
 +
 +Создадим объект этого класса uiExperiment: 
 +
 +  uiExperiment task;
 +
 +Будем считать, что пластинка бесконечна вдоль направлений x и y, и конечна вдоль направления z.
 +
 +Размер вычислительного объема вдоль напрвления z будет равным 3. Шаг по сетке равен 0.05. Мы используем поглощающие граничные условия вдоль z.
 +
 +  task.SetInternalSpace(Vector_3(0,0,0),Vector_3(0,0,3));
 +  task.SetResolution(0.05,1);
 +  task.SetBC(BC_PER,BC_PER,BC_PML);
 +
 +Мы указали нулевой размер вычислительного объема вдоль направлений x и y, поскольку наша задача фактически одномерна.
 +На самом деле в нашем случае толщина вдоль этих направлений будет равна двум шагам сетки по оставшемуся направлению z.
 +Граничные условия вдоль x и y будут периодическими, тем самым моделируя одномерность.
 +
 +Фактор Куранта должен быть не больше <tex>1/\sqrt{D}</tex>, где <tex>D</tex> это размерность задачи.
 +Мы используем фактор Куранта 1 (второй аргумент SetResolution), что допустимо в силу одномерности задачи.
 +Мы могли бы взять меньшее значение, но тогда шаг по времени был бы меньше.
 +
 +Световой импульс падает на пластинку снизу вдоль направления z.
 +Для генерации импульса зададим область Total Field в виде полупространства, ограниченное снизу плоскостью, которая перпендикулярна оси z и пересекает ее в точке 0.5.
 +
 +  task.AddTFSFPlane(INF,INF,0.5);
 +
 +Плоская волна распространяется вдоль оси z, вектор поляризации направлен вдоль оси x.
 +
 +  task.SetPlaneWave(Vector_3(0,0,1),Vector_3(1,0,0));
 +
 +В прошлом примере со сферой мы задавали явно временную зависимость падающего импульса в форме Беренгера с заданными параметрами.
 +В этот раз мы не будем ее задавать.
 +По умолчанию будет использоваться форма Беренгера с такими параметрами, чтобы спектр импульса хорошо покрывал длины волн, на которые приходится более 10 сеточных шагов.
 +
 +Поставим по одному детектору в область Scattered Field перед пластинкой и в область Total Field за пластинкой.
 +
 +  task.AddDetectorSet("h",Vector_3(0,0,0.25),Vector_3(0,0,2.75),iVector_3(1,1,2));
 +
 +Для расчета спектра потока вектора Пойнтинга отраженной и прошедшей волны мы поставим плоскости детекторов, перпендикулярные направлению z и пересекающие ось z в точках 0.25 и 2.75. 
 +
 +  task.AddRTASet("flux",2,0.25,2.75);
 +
 +Направления в функции AddRTASet задаются следующим образом: x – 0, y - 1, z – 2. Поэтому второй аргумент означает z.
 +
 +Поскольку поток вектора Пойнтинга падающей волны расчитывается программой автоматически по форме падающего импульса, для его расчета детекторов не требуется.
 +
 +Проведем расчет.
 +
 +  task.Calculate(100);
 +  task.Analyze();
 +
 +В результате граница TF/SF сгенерирует плоский импульс, который будет беспрепятственно распространяться вдоль направления z и в конце концов поглотится в PML. Детектор в зоне Total Field фиксирует распространяющийся импульс Беренгера
 +
 +  gnuplot> p 'h_r0001.d' u 1:5 w l
 +  
 +{{tutorial_plate:berenger_plate.png}}  
 +
 +Детектор в области Scattered Field фиксирует сигнал, отразившися от PML, который имеет амплитуду порядка 1е-4 от исходного сигнала (PML работает неидеально и что-то отражает). Однако это амлитуда мала и практически не влияет на результаты.
 +
 +  gnuplot> p 'h_r0000.d' u 1:5 w l
 +
 +{{tutorial_plate:berenger_refl.png}}  
 +
 +В файле flux.d будут представлены потоки вектора Пойнтинга отраженной, падающей и прошедшей волны соответственно. Их можно использовать для простроения спектра отражения и прохождения
 +
 +  gnuplot> p [][-.1:1.1] 'flux.d' u 1:($2/$3) w l t 'R', 'flux.d' u 1:($4/$3) w l t 'T'
 +  
 +{{tutorial_plate:rt_vac.png}}    
 +
 +Как мы видим, отражение равно нулю, а прохождение – единице.
 +
 +Зададим плоскопараллельную пластинку
 +
 +  emMedium M(2.25);
 +  M.AddPolarizability(emLorentz(1.1*(2*M_PI), 1e-5*(2*M_PI), .5));
 +  M.AddPolarizability(emLorentz(.5*(2*M_PI), .1*(2*M_PI), 2*1e-5));
 +
 +  task.AddObject(M,GetPlate(Vector_3(0,0,1),Vector_3(0,0,1),1));
 +
 +Пластинка задается функцией GetPlate, первый аргумент которой задает направление вдоль пластинки (направление z), второй точку, через которую проходит одна из сторон пластинки, а третий - толщину пластинки (она равна 1).
 +Пластинка сделана из дисперсной среды, описываемой двумя членами Лоренца.
 +
 +Повторив расчет, мы получим спектр отражения и прохождения для пластинки.
 +
 +  gnuplot> p 'flux.d' u 1:($2/$3) w l t 'R', 'flux.d' u 1:($4/$3) w l t 'T'
 +
 +{{tutorial_plate:rt_plate.png}}
 +
 +Мы видим, что результаты для маленьких частот неточны (прохождение больше 1). Эту ситуацию можно исправить, выбрав более широкий по времени импульс Беренгера, который покроет маленькие частоты.
 +
 +  task.SetSignal(new Berenger(1,2.5,.5,12.5));
 +  
 +Однако, в этом случае ухудшатся результаты для более высоких частот, в чем вы можете сами убедиться.
 +
 +Поглощение может быть получено как 1 - R - T.
 +
 +Расстояние между телами, границей TF/SF, PML и детекторами должно быть не менее 2 сеточных шагов. В нашем примере мы придерживались большего расстояния (5-20 сеточных шагов), поскольку иногда это может улучшить результат.
 +
 +=====Кремниевая подложка=====
 +
 +Ниже мы рассчитываем отражение от полубесконечной кремниевой подложки в случае нормального падения. Мы используем уже имеющийся код из предыдущей подсекции (задание граничных условий, области Total Field и т. д.)
 +
 +Для задания полубесконечной подложки мы можем воспользоваться функцией GetHalfSpace. В отличие от пластинки, у подложки нет второй замыкающей поверхности, и она бесконечно продолжается вдоль оси z. В модели FDTD она погружается в PML, что моделирует ее бесконечность.
 +
 +Материалом подложки будет кремний. Для получения соответствующего объекта emMedium можно воспользоваться функцией getSi, у которой имеется аргумент, по умолчанию равный единице. Этот аргумент задает, чему равна единица длины FDTD в микронах. Это необходимо для правильной калибровки коэффициентов, участвующих в расчете диэлектрической проницаемости (см. [[tutorial?&#sredy|как задавать среды в EMTL]]). В нашем случае единица FDTD будет равна 0.1 микрона.
 +
 +  task.AddObject(getSi(.1),GetHalfSpace(Vector_3(0,0,1), Vector_3(0,0,1)));
 +  
 +Параметры, используемые для апроксимации диэлектрической проницаемости кремния функцией getSi, были получены с помощью программы на MatLab, которую вы можете найти в разделе [[fitting]].
 +
 +По умолчанию в качестве PML используется тип UPML (Uniaxial PML). Однако, в него можно погружать только диэлектрические среды. Поскольку мы погружаем в PML дисперсную среду (кремний), мы должны использовать CPML (Convolutional PML). Он работает несколько медленнее, чем UPML, но является более универсальным.
 +
 +  task.SetPMLType(CPML_TYPE);  
 +
 +Построим зависимость отражения от кремниевой подложки от длины волны в мкм.
 +
 +  gnuplot> set xlabel 'wavelength, um'
 +  gnuplot> p [.3:1][0:1] 'flux.d' u (.1/$1):($2/$3) w l t 'R'    
 +
 +{{tutorial_plate:r_si.png}}
 +
 +Как видно, кремний достаточно сильно отражает в оптическом диапазоне. Для уменьшения отражения могут использоваться [[ar|антиотражающие текстурированные покрытия]].
 +=====Наклонное падение=====
 +
 +Итерационный метод для моделирования наклонного падения на периодическую структуру подразумевает проведение нескольких численных экспериментов FDTD, которые мы называем итерациями.
 +В течение каждой итерации происходит запись значений полей на периодических границах.
 +Это дает ключ к решению проблемы получения полей из будущего: поля в будущем известны из прошлой итерации, поскольку они были уже записаны.
 +При этом поля из прошлого могут браться из текущей итерации.
 +Для управления итерационным процессорм используется метод Total Field/Scattered Field (TF/SF).
 +Спустя некоторое число итераций граница TF/SF работает как применение периодических граничных условий с нужным сдвигом во времени.
 +
 +Наш метод описан в статье
 +[[http://www.opticsinfobase.org/abstract.cfm?uri=ol-33-13-1491|http]]{{:valuev_-_iterative_technique_for_analysis_of_periodic_structures_at_oblique_incidence_in_the_finite-difference_time-domain_method.pdf|PDF}}
 +
 +{{tutorial_periodic:geometry_obl.png?350}}\\
 +
 +Проиллюстрируем работу итерационного метода на примере расчета наклонного падения на фотонно-кристаллическую пластнику из упакованных в квадратную решетку сфер.
 +Мы используем метод Total Field/Scattered Field (TF/SF) для генерации плоской волны внутри области Total Field, которая заштрихована на рисунке. 
 +TF/SF граничные условия на границе 3 задают наклонно падающую волну, которая известна заранее. 
 +Идея итерационного метода заключается в применении дополнительной TF/SF генерации на границах <tex>{\bf x}_1</tex> и <tex>{\bf x}_2</tex>, используя записанные значения полей из точек-образов
 +<tex>{\bf x}_{1'}={\bf x}_1+{\bf a}</tex> и <tex>{\bf x}_{2'}={\bf x}_2-{\bf a}</tex>. 
 +Для отрицательного временного сдвига (поля из прошлого) можно брать поля из текущей итерации <tex>i</tex>, тогда как для положительного временного сдвига (поля из будещего) мы используем поля, записанные на предыдушей итерации:
 +
 +<tex>
 +{\bf F}_i({\bf x}_2,t)={\bf F}_i({\bf x}_{2'},t-\Delta t), \\
 +
 +{\bf F}_i({\bf x}_1,t)={\bf F}_{i-1}({\bf x}_{1'},t+\Delta t).
 +</tex>
 +
 +Расстояние между границами 1 и 2 больше периода <tex>a</tex> на некоторую величину <tex>\Delta x</tex> равную нескольким сеточным шагам. Это нужно для разделения границ TF/SF 1, 2 и точек-образов. 
 +
 +На первой итерации "падающая волна" на границах TF/SF 1, 2 есть нуль.
 +
 +Более подробно про итерационный метод можно прочитать здесь
 +[[http://www.opticsinfobase.org/abstract.cfm?uri=ol-33-13-1491|http]]{{:valuev_-_iterative_technique_for_analysis_of_periodic_structures_at_oblique_incidence_in_the_finite-difference_time-domain_method.pdf|PDF}}
 +
 +Для включения итерационного метода нужно использовать функцию SetOblique.
 +Первый параметр этой функции есть число итераций.
 +Второй параметр (назовем его n) отвечает за работу детекторов: детекторы работают и анализируется на каждой n-ой и на последней итерации.
 +Текстовые файлы детекторов с промежуточных итераций помечаются суффиксом  "_i" (i это номер итерации) после имени детектора.
 +
 +Ниже мы приводим код для моделирования наклонного падения в вакууме, рассматривается двумерный случай (размер вычислительного объема в направлении x есть нуль).
 +
 +<code cpp>
 +#include "uiexp.h"
 +
 +int main(int argc,char **argv){
 +
 +  emInit(argc,argv);
 +  uiExperiment task;
 +
 +  task.SetInternalSpace(Vector_3(.5,0,0),Vector_3(.5,2,3));
 +  task.SetResolution(.05);
 +  task.SetBC(BC_PER,BC_PER,BC_PML);
 +
 +  theta=45*M_PI/180.; // angle of incidence, radians. theta=0 corresponds to normal incidence
 +  Vector_3 k=Vector_3(0,sin(theta),cos(theta)); // light propagation direction
 +  Vector_3 E(1,0,0); // s polarization
 +//  Vector_3 E(0,cos(theta),-sin(theta)); // p polarization
 +
 +  task.AddTFSFPlane(INF,INF,.5);
 +  task.SetPlaneWave(k,E);
 +  task.SetOblique(10,1); // comment it for normal incidence
 +
 +  task.BuildDimensions();
 +
 +  Vector_3 p1,p2;
 +  task.GetSpace(p1,p2);
 +  Vector_3 mv1=p1+.05, mv2=p2-.05; // one mesh step from the calculated volume border
 +  mv1[0]=mv2[0]=(p1[0]+p2[0])/2;
 +  task.AddDetectorSet("movie",mv1,mv2,INT_INFTY,DET_TSTEP);
 +
 +  task.CalculateN(300);
 +
 +  task.Analyze();
 +
 +  return 0;
 +}
 +
 +</code>
 +
 +Заметим, что фактический размер вычислительного объема больше чем тот, который задан в аргументах функции SetInternalSpace.
 +Этот размер включает в себя PML и пространство слева от границы 1 и справа от границы 2.
 +Фактический размер вычислительного объема можно получить с помощью функции GetSize.
 +Для иллюстрации сходимости метода, мы поставили детекторы "movie" по всему вычислительному объему.
 +В течение численного эксперимента мы получим множество файлов moviexx_tyyyy.d, где xx это номер итерации, а yyyy это временной шаг.
 +Эти файлы можно использовать для построения эволюции распространения поля на разных итерациях.
 +
 +  gnuplot> set pm3d interpolate 2,2
 +  gnuplot> cbr=1.3
 +  gnuplot> set palette defined (-cbr 'blue', 0 'white', cbr 'red')
 +  gnuplot> set cbrange [-cbr:cbr]
 +  gnuplot> set view 0,0
 +  gnuplot> sp 'movie01_t0160.d' u 3:4:5 w pm3d
 +
 +{{tutorial_periodic:iter_obl.png?750}}\\
 +
 +Как можно видеть из рисунка, в нашем случае распределение поля на 160 временном шаге сходится к решению (плоский импульс в вакууме) за первые 5 итераций.
 +Для получения решения на более позднем временном шаге требуется больше итераций.
 +
 +=====Фотонно-кристаллическая пластинка=====
 +
 +Фотонно-кристаллическая пластинка это структура конечная вдоль одного направления и периодическая вдоль остальных направлений.
 +Ниже мы моделируем наклонное падение на такую пластинку состояющую из металлических (<tex>\epsilon=4</tex>, <tex>\sigma=2</tex>) сфер упакованных в квадратную решетку. 
 +Период решетки <tex>a=1</tex>, радиус сфер <tex>R=0.375</tex>, сеточный шаг равен 0.025.
 +
 +<code cpp>
 +#include "uiexp.h"
 +
 +int main(int argc,char **argv){
 +
 +  emInit(argc,argv);
 +
 +  uiExperiment task;
 +  task.SetInternalSpace(0,Vector_3(1,1,3));
 +  task.SetResolution(.025);
 +  task.SetBC(BC_PER,BC_PER,BC_PML);
 +
 +  valtype theta=45*M_PI/180.; // angle of incidence, radians. theta=0 corresponds to normal incidence
 +  Vector_3 k(0,sin(theta),cos(theta));
 +  Vector_3 E(1,0,0); // s polarization
 +//  Vector_3 E(0,cos(theta),-sin(theta)); // p polarization
 +
 +  task.AddTFSFPlane(INF,INF,.5);
 +  task.SetPlaneWave(k,E);
 +  task.SetOblique(100,1);
 +
 +  task.AddObject(emMedium(4,2),new Sphere(.375,Vector_3(.5,.5,1.5)));
 +
 +  task.AddDetectorSet("h",Vector_3(.5,.5,0.25),Vector_3(.5,.5,2.75),iVector_3(1,1,2));
 +  task.AddRTASet("flux",2,0.25,2.75);
 +
 +  task.CalculateN(2000);
 +
 +  task.Analyze();
 +
 +  return 0;
 +}
 +</code>
 +
 +Посмтроим геометрию в gnuplotе:
 +
 +  gnuplot> sp 'cell.pol' w l, 'med.pol' w l, 'pmlo.pol' w l, 'pmli.pol' w l, \
 +  gnuplot> 'tf.pol' w l, 'sg.vec' w vec, 'itf.pol' w l, 'unit.pol' w l
 +
 +{{tutorial_periodic:obl_gnu.png?350}}\\
 +
 +Согласно картинке выше, 'itf.pol' соответствует границам 1 и 2, а 'unit.pol' соответствует периодической ячейке ширины a (она ограничена границей между 1 and 2' и границей между 1' и 2).
 +
 +Сравним результаты для спектра прохождения, полученные с помощью FDTD и полуаналитического метода layer Korringa-Kohn-Rostoker (LKKR).
 +
 +{{tutorial_periodic:obl_pc_trans.png?450}}\\
 +
 +Как можно видеть из сравнения, результаты для больших частот сходятся быстрее.
 +При этом 5-10 итераций достаточно для получения спектра прохождения с нормальной точностью.
 +Для окончательной сходимости метода в нашем случае требуется порядка 100 итераций.
 +Вообще, число итераций, необходимое для сходимости метода, пропорционально числу временных шагов в численном эксперименте (в нашем примере это 2000).
 +
 +======Радиационное затухание источника======
 +
  
 
/var/www/fdtd.kintechlab.com/docs/data/pages/ru/tutorial.txt · Последние изменения: 2014/01/11 05:04 — deinega     Наверх