Моделирование диаграмм электромагнитного каротажа на графических процессорах

Вычислительные технологии

Том 13, № 6, 2008

Моделирование диаграмм электромагнитного каротажа на

графических процессорах

В.Н. Глинских, М.И. Эпов, И. Б. Ллвутин Институт нефтегазовой геологии и геофизики СО РАН, Новосибирск, Россия

e-mail: GlinskikhVN@ipgg.nsc.ru

The work is intended for development of a software-algorithmic base for expresssimulation and inversion of high-frequency electromagnetic logging data. Based on linear notions regarding characteristics of the electromagnetic field, the effective algorithm for simulation using graphics processors of personal computers (video cards) has been developed. Comparative assessments of processing power of logs for 2D geoelec-trical models of water saturated reservoirs for the central processing unit (CPU) and graphics processing unit (GPU) have been obtained.

Введение

При изучении геологической среды вокруг нефтегазовых скважин большую роль играет реконструкция пространственного распределения удельной электропроводности горных пород. Один из известных способов решения этой задачи состоит в регистрации откликов от вихревых токов в среде, возбуждаемых индуктивным источником в скважине. Их наиболее распространенными модификациями являются низкочастотный (20-250 кГц) индукционный каротаж с измерением реальной и мнимой частей ЭДС [1-3]. В последние десятилетия интенсивно развивается высокочастотный (0.8-14 МГц) электромагнитный каротаж (ВИКИЗ) с измерением относительных характеристик (разности фаз и затухания амплитуды) между сближенными приемными катушками [4-7].

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

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2008.

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

В отличие от центрального процессора (CPU) современные графические процессоры (GPU), специально спроектированные для очень быстрой обработки большого объема графических данных, способны исполнять параллельно тысячи потоков. В частности, графический процессор Nvidia GeForce 8600 GTS имеет теоретическую (пиковую) производительность 139 Gflops (число операций с плавающей точкой в 1 с) для чисел одинарной точности. Для сравнения, максимально возможная производительность процессора Intel Core 2 Duo 3 GHz — всего лишь 24 Gflops для чисел двойной точности и 48 Gflops — для одинарной точности. Следует отметить, что и полоса пропускания памяти графических процессоров гораздо выше, чем CPU. Например, пропускная способность памяти видеокарты Nvidia GeForce G280 составляет 141.7 Гбайта/с.

Использование графических процессоров для ускорения вычислений привлекает все больше внимание исследователей при решении современных задач в разных научных областях: от IT-технологий до медицины. Этому способствует специально разработанная технология Nvidia CUDA (Compute Unified Device Architecture), предназначенная для создания программного обеспечения для графических карт Nvidia [7]. Эта технология дает доступ к вычислительным ресурсам графического процессора через программный интерфейс для разрабатываемых приложений с большим объемом вычислений. При этом разработка программного обеспечения с помощью технологии Nvidia CUDA сильно упрощена благодаря использованию языка С с небольшими расширениями.

Все это указывает на то, что можно достичь ускорения расчетов при решении задач современной геоэлектрики с помощью вычислений на графических картах. В работе получены оценки производительности при расчетах диаграмм высокочастотного электромагнитного каротажа в типичных двумерных геоэлектрических моделях водонеф-тенасыщенных коллекторов при идентичных вычислениях на центральном процессоре Intel Core 2 Quad 2.4 GHz и видеокартах Nvidia GeForce 8600 GTS и G280.

1. Решение задачи

Рассмотрим двумерную геоэлектрическую осесимметричную модель, включающую пачку слоев с плоскопараллельными горизонтальными границами, пересеченную вертикальной скважиной кругового сечения. Каждый из пластов может содержать зону проникновения и/или окаймляющую зону, которые отделены друг от друга, от скважины и пласта коаксиально-цилиндрическими границами (рис. 1).

Каждая из областей, выделенная горизонтальными и вертикальными границами, характеризуется своими значениями электропроводности ají, j =1,...,Nz, l = 1,... , Nr (j), где Nz — число пластов, Nr (j) — число цилиндрических зон в соответствующем пласте. Будем полагать, что среда немагнитная ^jí = = 4п10-7 Гн/м. В остальной части среды электропроводность является функцией двух пространственных координат a(r,z). Предполагается, что электропроводность бурового раствора постоянная и равна ac.

Для решения двумерной прямой задачи высокочастотных электромагнитных зондирований в проводящих средах используем подход, изложенный в [8]. Здесь приведены полная линеаризованная постановка краевой задачи и ее решение — определение компонент электрического и магнитного полей на оси симметрии двумерной модели,

h

Рис. 1. Двумерная геоэлектрическая модель

с

г

с

11

z

1

22

23

г

г

21

22

z

2

31

32

г

31

разделенной на области с разной электропроводностью системой плоскопараллельных и коаксиально-цилиндрических границ. Этот подход состоит в представлении двумерного распределения электропроводности и электрического поля, возбуждаемого индуктивным источником, в виде суммы a(r, z) = ab(p) + z) и E = Eb + 5E. Здесь ab(p) и Eb — одномерное распределение электропроводности и электрическое поле в фоновой модели, а $a(r, z) и 5E — их относительно малые возмущения. Отметим, что ab зависит только от одной из пространственных координат (ab(r), ^b(z)) или является постоянной ab = ao (однородная среда). С учетом этого решение задачи определения электрического поля E(r0,z0|r,z) в двумерной проводящей среде сводится к решению следующего интегрального уравнения:

E

s

Здесь Eb(ro, zo|r, z) — электрическое поле в фоновой модели с распределением электропроводности ab(p), GE(r', z'|r, z) — функция Грина для электрического поля, ifc2(r', z') — возмущение квадрата волнового числа, связанное с относительно небольшими пространственными вариациями электропроводности, (ro,zo), (r, z) и (r',z') — координаты генераторной и приемной катушек и точки интегрирования. Осевая симметрия приводит к двумерному интегрированию по площади S неоднородности в плоскости ^ = const.

Магнитное поле и соответствующая функция Грина определяются исходя из выражения (1), а также H = --rotE и GH = --rotGE, и являются решением следующего интегрального уравнения:

H

s

Предположим, что пространственные вариации электропроводности приводят к изменению только фазы электрического или магнитного полей:

Е(го,гоМ) ~ Еь(го,гоМ)е^Е(г0'"0М), Н(го, *)|г, г) « Нь(го, ¿о|г, г)е^Н(го'"о|г'г).

Фазовые функции ^(го, го|г, г) строятся на основе решения тех же интегральных уравнений (1) и (2), при этом в подынтегральных выражениях полные поля заменяются на фоновые.

Отметим, что в силу симметрии среды и источника, расположенного на оси скважины, векторы электрического и магнитного полей содержат следующие ненулевые компоненты: Е = {0,Е^, 0} и Н = {Нг, 0, Н}. То же самое можно записать для соответствующих функций Грина: СЕ = {0, 0^, 0} и Он = {0^, 0,0^}. Тогда линейные представления для тангенциальной компоненты электрического поля и вертикальной компоненты магнитного поля записываются следующим образом:

/ —2п/ 5к2(т', г')0Е(т', г'|т, г)Е%(то, го|т', г')^Б

Е^(го,го|г, г) к Е^(то,го|т, г)ехр

v v ' i '

"

\

Е£

Н(то,го|т, г) = 7777--т^тЕ^то, го|т, г).

1 1 д т дт

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

Линеаризуем разность фаз ЭДС между двумя сближенными приемными катушками (Д<£ = — Поскольку наведенная в катушке ЭДС пропорциональна вертикальной компоненте магнитного поля е = , то

ТтН

= агс1§^еи:, г =1,2'

Здесь и далее индекс г для вертикальной компоненты магнитного поля опущен. Так как расстояние между сближенными катушками ДЬ много меньше их расстояния до генераторной катушки Ь (т.е. ^ 1), то разность фаз можно представить через вер-

дю

тикальную производную фазы Дю к —— ДЬ. Согласно (4) выражение для производной

дг

фазы запишется следующим образом

дЮ = 1 ^ТтИд РеИд 1тН ^

где А = у^Ке2Н + Тт2Н — амплитуда магнитного поля.

Линейное приближение для вертикальной компоненты магнитного поля на основе (3) в общем виде можно записать

( 2п / 8НйБ^

Н к Ньехр

5к2

\

Нь

где 8Н — подынтегральная функция двумерного интеграла по области возмущения электропроводности, описываемая в (3). Основываясь на этом выражении, фазу вертикальной компоненты магнитного поля Н запишем в следующем виде:

( 2п / 8НйБ^

Ю к юь + Тт

Як-2 Ь

дк НЬ

\

Дифференцируя это выражение по г в соответствии с (5) и выделяя во втором слагаемом мнимую часть, получаем линейное представление разности фаз:

Д^ « Д^ь + + /2)ДЬ.

— разность фаз в фоновой среде, — возмущения электропроводности, а ве-и /2, являющиеся двумерными интегралами по определяются следующим

/! = 2п / (ве«^ ^ + ^ д

г / д 1шЯь 1шЯь 5 А

/2 = 2^ (1шгнддт ^ + ^ а; мяЫ

5 у 7

Отметим, что точность линейного представления зависит от выбора фоновой модели, размеров неоднородности и относительных контрастов электропроводностей в среде. Сравнительный анализ численного моделирования линейных аппроксимаций электромагнитного поля и его относительных характеристик, оценки точности, эффективных приемов улучшения, выбора фоновой модели и иных процедур проведен в [8]. В качестве фонового используем электромагнитное поле горизонтального витка с единичным током в двухслойной цилиндрической среде [9].

Здесь Д^ь личины Д образом:

2. Схема интегрирования

Для построения эффективной схемы интегрирования и расчета интегралов /! и /2 детально проанализированы их подынтегральные функции (7). Эти коэффициенты-функции под интегралами по области возмущения электропроводности представлены вещественными полиномами по переменной кь и зависят от длины зонда и координат точки интегрирования.

'¡1+1

пп

• • •

г

п

г

г

п

в

п

г

а

А

Рис. 2. Схема интегрирования

Организуем интегрирование коэффициентов-функций и вычисление интегралов Д и 12 в соответствии с параметризацией модели по областям S возмущения электропроводности 8uji в порядке j = 1,... , Nz, l = 1,... , Nr (j). Общая схема интегрирования приведена на рис. 2. Обратим внимание на то, что области 8oji могут иметь существенно разные размеры, поэтому для уверенной сходимости интегрирования необходимо унифицировать сетку интегрирования. Для этого каждая область 8aji дополнительно разбивается по z и r на области с характерным размером 8, определяемым величиной скин-слоя в фоновой модели. В этом случае каждая область возмущения электропроводности характеризуется своим числом разбиений nz и nr. Для интегрирования в подобласти 8 используем метод Гаусса с ng-точечной квадратурой. Расчеты проведем с использованием 3-, 6-, 9- и 16-точечных квадратур Гаусса (ng = 3, 6, 9, 16).

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

3. Вычислительный алгоритм на графической карте

Обратим внимание, что задача линеаризации относительных характеристик электромагнитного поля, приводящая к вещественным величинам, является самостоятельной и представляет особый интерес при решении прямых и обратных задач. Отметим, что расчеты, в частности на графической карте Nvidia GeForce 8600 GTS, сделаны для вещественных чисел. Карта базируется на чипе мультипроцессора G84 и на ней возможны вычисления с плавающей точкой одинарной точности.

Доступ к вычислительным ресурсам графического процессора осуществляется через специально разработанную для этих целей технологию Nvidia CUDA. Предназначенная для разработки программного обеспечения и алгоритмических приложений для графических карт, она дает возможность доступа к набору инструкций графического ускорителя, средствам управления его памятью и позволяет организовать на нем сложные параллельные вычисления.

К основным особенностям архитектуры CUDA можно отнести следующие. Графический процессор рассматривается как устройство, способное параллельно исполнять множество потоков. При разработке программного обеспечения (реализации алгоритма) в приложении выделяется главная его часть, называемая ядром. Ядро — часть приложения, параллельно исполняемая потоками непосредственно на графическом процессоре. Вычислительные потоки организованы в блоки. Объединенные в блок потоки взаимодействуют с помощью разделяемой памяти. Расположенная непосредственно на чипе мультипроцессора графической карты разделяемая память имеет небольшой объем (16 Кбайт для карт на ядре G84) и эффективна по быстродействию при операциях чтения и записи, обладая высокой пропускной способностью. Указанные блоки потоков имеют свою размерность. При этом блоки одной размерности объединяются в грид и исполняют одно и то же ядро. Каждый поток и блок имеют свои индексы, являющиеся номером потока внутри блока и номером блока внутри грида соответственно. Поэтому как блоки, так и гриды могут рассматриваться как массивы потоков и блоков.

blockID.r -

(2,1) (2,2) (2,3) adID.r -

(3,1) (3,2) N DI J 1 thre

(4,1) _

Рис. 3. Модель алгоритма на CUDA

Существует возможность сконфигурировать блок как одно-, дву- или трехмерный массив, а грид — как одно- или двумерный.

С помощью технологии Nvidia CUDA алгоритм вычисления значений интегралов 1\ и 12 по приведенной выше схеме интегрирования естественным образом ложится на программную модель CUDA и не требует специфической организации памяти для хранения данных и промежуточных расчетов, что делает его простым в понимании и реализации.

В разработанном алгоритме для вычисления значений интегралов Д и 12 область интегрирования разбивается на подобласти общим числом nz х nr, в каждой из которых методом Гаусса вычисляется двумерный интеграл. Для ng-точечной квадратуры Гаусса в каждой подобласти необходимо вычислить и просуммировать ng х ng чисел. Соответственно для вычисления интеграла необходимо просуммировать его значения в каждой подобласти. Таким образом, ядро программы параметризуется числом подобластей и числом точек квадратуры Гаусса, поэтому в алгоритме используются двумерные блоки и гриды, а соответственно и двухкомпонентные номера блоков (blockID.z и blockID.r) и потоков (threadlD.z и threadlD.r). Каждая подобласть интегрирования обрабатывается одним блоком, в свою очередь каждый поток вычисляет значение подынтегральной функции в одной точке квадратуры Гаусса. Результаты вычислений каждого потока хранятся в разделяемой памяти. Конфигурация ядра программы вычислений интеграла отражает разбиение области интегрирования (рис. 3).

Обратим внимание на вычислительные возможности используемых графических карт. Так, графическая карта Nvidia GeForce 8600 GTS позволяет организовать максимальное число потоков в блоке, равное 512, при этом максимальный объем разделяемой памяти на блок составляет 4096 байт. В разработанном алгоритме при реализации интегрирования в подобласти с использованием 16-точечной квадратуры Гаусса в блоке организуется 256 потоков и соответственно для хранения результатов вычислений в блоке используется 1024 байта разделяемой памяти (4 байта для числа с плавающей точкой одинарной точности).

Сравнительный анализ диаграмм

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

для терригенных разрезов нефтегазовых месторождений Западной Сибири, использовав значения кажущейся удельной электропроводности для зондов ВИКИЗ: 1К05 — длина 0.5 м, частота 14 МГц; 1К07 — 0.7 м, 7 МГц; 1К10 — 1.0 м, 3.5 МГц; 1К14 — 1.4 м, 1.75 МГц; 1К20 — 2.0 м, 875 кГц.

Тестирование алгоритмов и определение пределов применимости приближенного подхода к моделированию диаграмм на основе выражения (6) проведено на модели водоплавающего нефтенасыщенного коллектора. Она включает скважину (ас =1/2 См/м, гс = 0.1 м), нефтенасыщенный коллектор с зоной проникновения и окаймляющей зоной (к = 2 м, а23 =1/8 См/м, г22 = 0.6 м, а22 = 1/4 См/м, г21 = 0.4 м, а21 = 1/20 См/м), подстилающий водонасыщенный коллектор с зоной проникновения (а32 = 1/5 См/м, г31 = 0.5 м, а31 = 1/25 См/м) и перекрывающие глинистые отложения (ац = 1/4 См/м) (см. рис. 1). На рис. 4 приведены синтетические диаграммы кажущейся электропроводности ак для зондов ВИКИЗ. Для верификации расчетов была использована программа вычисления диаграмм на основе конечно-разностного алгоритма [10]. Как видно, результаты приближенного численного моделирования хорошо согласуются с конечно-разностным методом. На всем интервале при 4-кратном относительном контрасте относительная погрешность для зондов не превышает 5 %. Так, на интервалах нефте-и водонасыщенных коллекторов относительная погрешность для зондов 1К05 и 1К07 составляет около 3 и 2.5 %, для 1К10.. .1К20 — 1-2 %. На интервале глинистой покрышки значения относительной погрешности несколько увеличиваются для всех зондов. Максимальные относительные погрешности связаны с короткими высокочастотными зондами, поскольку линейное приближение дает приемлемую точность при малых относительных контрастах и сравнительно низких частотах. Для длинных зондов относительные погрешности не превышают 5 % даже при больших относительных контрастах. Аппроксимация (3) как способ коррекции фазовой функции позволяет получить более точные значения кажущихся электропроводностей. Предположение о том, что пространственные вариации электропроводности в высокочастотном диапазоне приводят к изменению только фазы компонент электрического или магнитного полей, вполне

1 2 Глубина, м

Рис. 4. Синтетические диаграммы в модели водоплавающего нефтенасыщенного коллектора, перекрытого глинами: 1 — приближенный, 2 — конечно-разностный подходы

справедливо. Выполненный большой объем численного моделирования и сравнительный анализ диаграмм в двумерных моделях терригенных коллекторов показал высокую эффективность предложенного подхода. Применительно к высокочастотному электромагнитному каротажу использование линейного представления разности фаз (6) позволяет проводить вычисления с относительной погрешностью не более 5 % в достаточно широком диапазоне модельных параметров.

4. Сравнение производительности при вычислениях

Получены оценки производительности для высокочастотного электромагнитного каротажа в двумерной модели нефтенасыщенного коллектора при вычислениях на центральном процессоре Intel Core 2 Quad 2.4 GHz и графических картах Nvidia GeForce 8600 GTS и G280. Для оценки времени расчетов были использованы 3-, 6-, 9- и 16-точечные квадратуры Гаусса. При этом число подобластей nz х nr, определяемое на этапе оценки характерного масштаба 5 с учетом геометрических размеров области возмущения электропроводности Saji, было постоянным. Укажем, что алгоритмы, реализованные при вычислениях на центральном процессоре и графических картах, полностью идентичны. На рис. 5 показано время вычислений диаграмм кажущейся электропроводности а^ для пяти зондов ВИКИЗ по вертикальному профилю (61 точка) на интервале нефте-насыщенного коллектора при интегрировании на двумерной сетке с различным числом квадратур Гаусса. Как видно, при использовании GPU для вычислений удается достичь увеличения производительности в 10-30 раз по сравнению с идентичными расчетами на CPU. Так, время вычисления диаграмм при использовании 6-точечной квадратуры Гаусса составляет около 1 с на CPU и 0.5 с на видеокартах. Десятикратное увеличение производительности достигается при использовании 9-точечной квадратуры на GPU 8600, при этом на GPU 280 — чуть более 16 раз. При использовании 16-точечной квадратуры Гаусса время вычисления диаграмм на CPU составляет около 27 с, при этом на GPU 280 — менее 1 с. Обратим внимание на то, что существенно большее время

на CPU и GPU

30

CPU

25

0

5

GPU8600

GPU280

3

6

9

16

Число точек квадратуры Гаусса

Рис. 5. Производительность вычислений на CPU и GPU

требуется для расчета диаграмм c использованием конечно-разностных алгоритмов на CPU, которое составляет более 30 мин. Безусловно, такое прямое сравнение ресурсо-емкости расчетов некорректно, поскольку оно связано с влиянием числа узлов сетки интегрирования и разностной схемы. Можно ожидать существенный прирост производительности конечно-разностных алгоритмов при расчетах на графических процессорах, что открывает новые возможности для оперативного моделирования и инверсии данных каротажа.

Заключение

На основе специализированных вычислений на графических процессорах персональных компьютеров (видеокартах) разработан эффективный алгоритм моделирования диаграмм высокочастотных электромагнитных зондирований в нефтегазовых скважинах. Получены оценки производительности при расчетах диаграмм электромагнитного каротажа в двумерных геоэлектрических моделях при вычислениях на центральном процессоре и графических картах. Показано, что при использовании GPU для вычислений удается достичь увеличения производительности в 10-30 раз по сравнению с идентичными расчетами на CPU.

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

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

Список литературы

[1] Doll H.G. Introduction to induction logging and application to logging of wells drilled with oil base mud //J. Petroleum Technology. 1946. Vol. 6, N 4. P. 148-162.

[2] Кауфман А.А. Теория индукционного каротажа. Новосибирск: Наука, Сиб. отд-ние, 1965. 236 с.

[3] Плюснин М.И. Индукционный каротаж. М.: Недра, 1968. 142 c.

[4] Антонов Ю.Н. Изопараметрическое каротажное зондирование (обоснование ВИКИЗ) // Геология и геофизика. 1980. № 6. С. 81-91.

[5] Технология исследования нефтегазовых скважин на основе ВИКИЗ. Методическое руководство / Ред. М.И. Эпов, Ю.Н. Антонов. Новосибирск: НИЦ ОИГГМ СО РАН. Изд-во СО РАН. 2000. 121 с.

[6] Даев Д.С. Высокочастотные электромагнитные методы исследования скважин. М.: Недра. 1974. 191 с.

[7] CUDA 2.0 Programming Guide. http://developer.download.nvidia.com/compute/cuda/2_0/ docs/NVIDIA_CUDA_Programming_Guide_2.0.pdf

[8] Эпов М.И., Глинских В.Н. Электромагнитный каротаж: моделирование и инверсия. Новосибирск: Академическое изд-во "Гео", 2005. 98 с.

[9] ТАБАровский Л.А. Применение метода интегральных уравнений в задачах геоэлектрики. Новосибирск: Наука, Сиб. отд-ние, 1975. 140 с.

[10] Эпов М.И., Мартаков С.В. Прямые двумерные задачи электромагнитного каротажа // Геология и геофизика. 1999. № 2. С. 249-254.

Поступила в редакцию 17 ноября 2008 г.