RU (546769), страница 2

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 2)

Уравнение (26) может быть переписано в форме, похожей на пропагатор (23):

. (27)

Такая трансформация была названа «звездным произведением» [7]. Так как дифференциальное уравнение (23) для пропагатора известно, имеем возможность получить одномерную задачу с начальными условиями матричного уравнения Рикатти для матрицы элементов рассеивателей [7] например, отражение:

, (28)

и аналогично для передачи .

Легко видеть, что уравнение (28) эквивалентно уравнению Амбарцумяна-Чандрасекара [8], полученного по принципу инвариантных вложений. Таким образом, после удаления анизотропной части и дискретизации уравнения мы получим единственное аналитическое решение в матричной форме (26). В этом смысле метод последовательных порядков рассеяния является итерационным методом для нахождения обратной матрицы в (26), а метод Монте-Карло – стохастический.

Реализация предложенного решения (26) в качестве алгоритма включает вычисление следующих компонент: нулей и весов квадратурной формулы для дискретизации ВУПИ (18), функций источника (4), собственных векторов и значений матрицы системы (24), произведений матриц (26).

Практическая реализация указанного алгоритма зависит в основном от размеров матриц в (26). Заметим, что (26) получено для реальной части, таким образом, размер зависит в основном не от средних оптических параметров и граничных условий, а от устраненной мнимой части решения. В общем случае, . Однако, в случае правильного устранения мнимой части решения, когда гладкая реальная часть находится в окрестности изотропного углового распределения, вполне возможно, что . Это значительно повышает производительность алгоритма.

4. Ликвидация мнимой части решения

Проанализируем и сравним все основные известные методы ликвидации мнимой части решения с точки зрения уменьшения размеров M и N матрицы и повышения эффективности поиска решения (26).

4.1 Ликвидация прямых излучений без рассеивания

Впервые Эддингтон предложил устранить мнимую часть путем вычислений прямых не рассеивающихся излучений на основе закона Бугера, Милн развил эту теорию, которая подтвердилась в работах Чандрасекара [8]. Таким образом, выражение для прямого не рассеивающегося излучения на плоскости совпала с плоскостью падения, определяющейся по формуле:

(29)

где – вектор параметров Стокса падающего излучения, – нормаль к лучу излучения, p, q – линейная и круговая степени поляризации луча, – азимут поляризации.

В этом случае функция источника в уравнении для реальной части принимает вид:

(30)

В этих условиях излучение рассеивается на малые углы и незначительно отличается от прямого излучения, и этот метод становится неэффективным. Большое рассеяние приводит к значительному увеличению числа членов разложения матрицы рассеяния в ряд по обобщенным сферическим функциям K, что, соответственно, увеличивает размеры матрицы N и M в решении (26).

4.2. Дельта метод

Некоторые работы были выполнены для решения проблем мнимых решений. Были использованы различные алгоритмы, основанные на усечении функции рассеивания фаз. Дельта метод наиболее эффективный среди них; он может быть обобщен и в векторном случае [9]. Как и в скалярном случае, функция рассеивания представляется в виде суммы дельта-функций и гладкой части.

(31)

Здесь f является мнимой частью рассеивания; – единичная матрица.

Уравнение (31) для греческой матрицы (11) можно записать как

(32)

Этот этап матричных преобразований приводит к масштабным преобразованиям оптической глубины и коэффициента альбедо рассеивания.

. (33)

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

4.3. TMS и IMS-методы

Для устранения упомянутых выше проблем в решении ВУПИ Накаджима и Тенака [10] обратились к идее определения мнимой части в решении, основываясь на аппроксимации аналитического представления углового распределения параметров Стокса первого и второго порядков рассеивания:

(34)

Здесь первые два порядка удовлетворяют уравнениям:

(35)

(36)

Уравнения (35) и (36) решаются аналитически. Однако решение уравнение (36) вычисляется с помощью тройного интеграла, и время выполнения его расчета может превысить время решения исходной задачи (26). По этой причине, как правило, используется приближение, описанное разбиением параметров Стокса в окрестности инцидентного направления. Наиболее эффективным способом решения является приближение малыми углами, где дисперсия рассеиваемых и не рассеиваемых лучей не принимается во внимание. Это приближение эквивалентно изменению параметра в формулах (35) и (36) на . Наиболее известные коды, использующие как дельта метод, так и TMS метод, – DISORT [11] в скалярном случае и Pstar [12] – в векторном.

4.4. Метод малых угловых изменений сферических гармоник

Может быть показано, что естественной эволюцией TMS-метода является включение в всех порядков рассеивания в малом угловом приближении [6]. Этот подход описан в [6] на основе метода малых угловых модификаций сферических гармоник.

В этом случае мнимая часть углового распределения параметров Стокса в CP-представлении расширен сферическими функциями с инцидентным направлением :

(37)

здесь , – азимут в системе относительно .

Спектр сильного мнимого углового распределения – гладкая функция гармонических индексов k. По этой причине можем получить дифференциальное уравнение [6] для коэффициентов разложения

(38)

решением которого является

(39)

После преобразования SP-представления МСГ [6] решение может быть представлено в системе координат вокруг нормали к границе слоя:

(40)

где .

После подстановки выражения (40) в (4) и после некоторых преобразований [6], получаем выражение для функции источника

(41)

где , .

Такой подход приводит реальную часть решения почти мнимой функцией. Таким образом, M может быть небольшой, N практически не зависит от K. Этот метод ликвидации мнимой части решения реализован в двух кодах: МДО в скалярном и векторном случаях. Оба кода реализованы на языках Fortran и MATLAB [13].

5. Численная реализация

Прежде всего отметим, что численная реализация решения (26) основана на решении основных задач линейной алгебры. Существуют несколько библиотек, в которых собраны необходимые подпрограммы: LAPACK, IMSL, NAG и так далее. Math Kernel Library (MKL) представляет собой пакет, который специально оптимизирован для процессоров Intel (существуют версии MKL для Windows и для Linux систем). Отличное преимущество получают от использования поддержки многоядерных и многопроцессорных систем и автоматического распараллеливания. Это способно значительно увеличить скорость вычислений. Существует вектор реализаций элементарной операции (например, квадратный корень из элемента массива) в параллельном режиме. Специальные математические пакеты программ такие, как Matlab, Maple, Mathematica используют MKL, играя роль библиотечной оболочки.

Перспективный путь для ускорения кода – применение графических процессоров общего назначения (GPUs). Технология CUDA является вычислительным двигателем в графических процессорах NVIDIA. Она дает доступ к набору команд вычислительных элементов видеокарты. Стандартные процедуры LAPACK могут быть ускорены CUDA. Вместо того, чтобы перепрограммировать процедуры LAPACK, желательно применять специальные CUDA-библиотеки. Кроме того, MATLAB 2010b поддерживает nVIDIA CUDA-совместимые процессоры и для использования GPU-вычислительных возможностей не требует знаний и опыта в CUDA.

Значительное ускорение и экономия памяти могут быть достигнуты при работе с разреженными матрицами. Основная идея заключается в следующем: при работе с матрицей, содержащей много нулевых элементов, имеет смысл хранить только ненулевые элементы и дополнительную информацию, которая может быть использована для восстановления индексов ненулевых или нулевых элементов. Существует несколько форматов для хранения разреженных матриц: формат сжатия разреженных строк (CSR), формат сжатия разреженных столбцов (CSC) и координированный формат. Выбор формата зависит от типа операции, например, CSC формат удобно использовать в случаях, когда элементы разреженной матрицы расположены в каждом столбце. Использование разреженных матриц может быть эффективным для вычисления мнимой части (49) и модифицированной функции источника.

6. Результаты и обсуждение

Были проанализированы влияния различных аппаратно-программных средств для эффективной реализации решения (26) с ликвидацией мнимой части на основе MSH на примере кода MDOM и MVDOM. Сравнение времени выполнения для различных режимов компиляции и вычисления представлены в таблице 1 для двух тестов (тест 1: N=101, K=500, M=32; тест 2: N=101, K=1000, M=32). Тесты проводились на системе Ubuntu 10.04, Intel Core 2 Duo 3 ГГц 2 Гб ОЗУ, Intel Fortran Compiler 11.1 c MKL 10.2. Использовались два компилятора gfortran и ifort, использовалась технология оптимизации при работе с разреженными матрицами.

Обратим внимание, что MKL использует все доступные ядра вычислительной системы (в данном случае – 2 ядра) и предоставляет результат в половину времени. Технология разреженных матриц, примененная при вычислении мнимой части, значительно уменьшила время выполнения. В связи с разреженностью матриц, двумерные массивы сводятся к одномерным. Таким образом, время выполнения пропорционально K вместо .

Мы увеличили ускорение примерно на 20% за счет умножения матриц на nVIDIA GeForce 480 GTX GPU. Вычисления на графических процессорах дают преимущества только для массивов большой размерности, в противном случае более предпочтительными являются вычисления на центральном процессоре. Инструменты профайлера показывают, что вычисление собственных векторов и значений занимают половину времени выполнения. К сожалению, подпрограммы для этих задач не реализована в инструментах Matlab для графических процессоров. Размеры матриц значительно сокращены в связи с тем, что для ликвидации мнимого решения использовался МСГ. На практике N меньше 300. Однако вычисление поляризационного спектра длины волны также может быть реализовано средствами CUDA технологии. Возможно использовать преимущества CUDA для более чем одной длины волны.

Характеристики

Список файлов перевода