Содержание к диссертации
Введение
Глава 1. Численное моделирование солнечной конвекции и рассеяние звуковых волн 27
1.1 Численное моделирование двумерной сжимаемой солнечной конвекции и супергрануляции 27
1.2 Моделирование рассеяния акустических волн солнечным пятном 41
Глава 2. Применение регуляризирующих алгоритмов для анализа солнечных наблюдений 55
2.1 Восстановление угловой скорости внутреннего дифференциального вращения Солнца по гелиоссйсмологическим данным 55
2.2 Определение аппаратной функции телескопа и восстановление изображений солнечных протуберанцев 61
Глава 3. Численное моделирование и анализ экспериментальных данных в задачах экологии 71
3.1 Восстановление функции распределения частиц по радиусам из седиментационной кривой 71
3.2 Численное моделирование седиментации полидисперсной суспензии в присутствии двумерной сжимаемой конвекции 79
3.3 Восстановление скорости роста и функции плотности вероятности из зашумленных экспериментальных данных 93
Заключение 112
Библиографический список использованной литературы 116
- Моделирование рассеяния акустических волн солнечным пятном
- Определение аппаратной функции телескопа и восстановление изображений солнечных протуберанцев
- Численное моделирование седиментации полидисперсной суспензии в присутствии двумерной сжимаемой конвекции
- Восстановление скорости роста и функции плотности вероятности из зашумленных экспериментальных данных
Введение к работе
Любое научное исследование должно последовательно ответить на три вопроса: "Что?", "Как?" и "Почему?". Первый этап соответствует сбору информации. Второй этап связан с обработкой данных. Третий этап способствует выяснению причин, породивших исследуемое явление, и здесь не обойтись без построения моделей. В предлагаемой диссертации используется комплексный подход к научному исследованию, охватывающий этапы моделирования и обработки экспериментальных данных на примере решения различных задач. Так как автор работает в Солнечном отделе Крымской астрофизической обсерватории и занимается численным моделированием конвекции и обработкой наблюдений, то большинство задач, решаемых в диссертации, так или иначе, связаны с физикой Солнца, конвекцией, дифференциальным вращением и гелиосейсмо-логией. Все эти задачи тесно связаны друг с другом.
Актуальность темы. Дифференциальное вращение вызывается взаимодействием вращения и конвекции, что приводит к анизотропии рейнольдсовых напряжений. Солнечные пятна, по-видимому, так же являются результатом конвективных движений — нисходящих потоков холодного вещества, быстро остывающего в областях с сильным магнитным полем и представляют собой самоорганизующиеся магнито-конвективные структуры. Однако, динамика под-фотосферных структур и потоков вещества еще недостаточно хорошо изучена. Прямое наблюдение невозможно из-за непрозрачности вещества, поэтому для диагностики внутренних слоев Солнца приходится решать обратные задачи акустики, применяя технику время-расстояние, используемую в гелиоссисмо-логии. В целом, проблема солнечной конвекции и се влияние на вращение и структуру магнитного поля Солнца еще не решена. В диссертации развиваются отдельные модели и методы обработки наблюдательных данных для решения этой фундаментальной проблемы солнечной физики и астрофизики.
Подобные модели и методы могут быть применены для моделирования конвекции и интерпретации данных, полученных в лабораторных условиях. Рассмотрение этих задач позволяет более полно выявить характерные отличия и физические закономерности конвективных процессов, расширить область применения численных методов, и, кроме того, получить результаты важные для экологических приложений. В частности, в диссертации выявлено, что конвекция на Солнце и в лабораторных условиях носит качественно различный характер. Если в лабораторных условиях конвекция развивается в виде хорошо известных ячеек Бенара, в которых горячее вещество поднимается в центре ячеек и опускается на границах, то на Солнце структура конвекции определяется высокоскоростными структурами опускающегося вещества. Образно говоря, если в лабораторных условиях конвекция развивается «снизу вверх», то на Солнце — наоборот «сверху вниз».
Среди других интересных аналогий можно отметить процесс оседания пыли в среде с конвективными движениями при формировании солнечной системы и седиментации в лабораторных и промышленных условиях для разделения суспензий. Конвекция противодействует гравитационному оседанию, действуя по-разному на различные пылевые фракции и элементы. Это приводит к разделению среды по физическим свойствам. Большое внимание в современной физике Солнца уделяется изучению структуры и динамики «тахоклина» — переходной зоны между радиативным ядром и зоной турбулентной конвекции. Считается, что в этой зоне протекают процессы солнечного динамо. Гелио-ссйсмологические измерения указывают на существенное отклонение в значениях скорости звука в этой зоне по сравнению со стандартной моделью. Это может быть объяснено процессами гравитационного оседания гелия и тяжелых элементов. Пока эти процессы описываются с использованием феноменологических моделей турбулентной диффузии (Proffitt & Michaud, 1991; Christensen-Dalsgaard et al., 1993), однако, результаты диссертации по моделированию седиментации указывают, что на следующем этапе возможно применение прямого численного моделирования тахоклина.
Солнце — это единственная звезда достаточно близкая к нам, чтобы можно было различать детали конвективных движений. С одной стороны, Солнце является лабораторией для тестирования моделей звездной конвекции, с другой — исследование солнечной конвекции имеет самостоятельный интерес. Динамика конвективной области Солнца является определяющей для таких глобальных явлений, как дифференциальное вращение, солнечное динамо и 11-ти летний цикл солнечной активности, возбуждение 5-мин. осцилляции. На поверхности мы видим проявление конвекции, как минимум, на двух различных масштабах — грануляции и супергрануляции. В то время, как механизм возникновения гранул хорошо изучен — это приповерхностная турбулентная конвекция, супергрануляция до сих пор, начиная с ранних наблюдений (Hart, 1954), остается загадкой. Харт (Hart, 1956) описал пространственные структуры с характерным размером 26 Мм и средними скоростями порядка 0.3 км/с, однако отклонил гипотезу о конвективной природе этих образований из-за слишком болыпого масштаба. Симон и Лейтон (Simon & Leighton, 1964) сорок лет назад сформулировали гипотезу, которой придерживаются и современные исследователи, что за возникновение супер грануляции ответственны зоны ионизации гелия.
Грануляционная картина простирается вглубь не более, чем на 2+3 Мм. Однако, глубина супергрануляционного слоя в настоящее время не известна. Существуют предположения, что в сильно стратифицированной атмосфере вертикальный размер супергранул может намного превосходить их горизонтальный размер (Parker, 1973). Вычисляя коэффициенты корреляции между горизонтальной и вертикальной компонентами скорости вещества на различных глубинах, Дюваль (Duvall, 1998) получил оценку глубины супергрануляционного слоя в 8 Мм. Жао и Косовичев (Zhao & Kosovichev, 2003) с помощью усовершенствованного алгоритма многоканальной деконволюции установили наличие сходящихся потоков в супергрануляционной ячейке на глубине около 10 Мм и оценили среднюю глубину супергрануляционного слоя в 15 Мм.
По понятным причинам, для упрощения аналитического или полуаналитического описания, супергрануляцию рассматривают в рамках ламинарной конвекции, не смотря на то, что движения вещества в атмосфере Солнца должны быть исключительно высокотурбулентными и существенно нестационарными, как отмечалось Симоном и Лейтоном (Simon & Leighton, 1964). Недавние результаты численного моделирования конвекции в стратифицированной среде при больших числах Рэлея выявили исключительно сложную структуру движений вещества. В настоящее время не вызывает сомнения, что перенос тепла и импульса в высокотурбулентной солнечной конвекции регулируется серией "плюмов" — узко локализованных высокоскоростных потоков относительно холодного вещества, движущегося вниз от поверхности (Julien et al., 1996; Brummel et al, 1996). Динамика отдельных "плюмов" определяется сильным вихревым взаимодействием с соседними "плюмами" (Julien et al., 1996). Раст и Плонер (Rast, 2002; Ploner et al., 2000) считают, что супергрануляционный масштаб возникает в результате взаимодействия и наложения отдельных грануляционных "плюмов". Похожая модель была предложена (Rieutord et al., 2000), в рамках которой супергрануляция является проявлением нелинейной крупномасштабной неустойчивости грануляционной картины, приводимой в действие расширяющимися гранулами. В обоих этих моделях размер супергранул не является характерным масштабом конвекции, а глубина супергрануляционного слоя определяется исключительно глубиной проникновения "плюмов", пока они остаются стабильными.
До сих пор не ясен вопрос, существуют ли промежуточные характерные масштабы: мезогранулы (November et al., 1981) и гигантские ячейки (Beck et al., 1998). В работе (Hathaway et al., 2000) авторы приходят к выводу, что мезогранулы отсутствуют как отдельный характерный масштаб солнечной конвекции. Ячеистая структура с характерным размером мезогранул присутствует, но объясняется вкладом либо больших гранул, либо малых супергранул. Теория длины пути перемешивания говорит, что для каждого фиксированного уровня существует только один характерный размер конвективных ячеек, связанный со шкалой по давлению на этом уровне. Существование дискретш>іх характерных масштабов дает указания на то, что солнечная конвекция может быть обусловлена чем-то большим, нежели простой стратификацией.
Рядом авторов отмечаются необычные свойства супергрануляции. Долготная кросс-корреляция супергрануляционной картины потоков свидетельствует о том, что она (картина) вращается быстрее плазмы на поверхности Солнца (Duvall, 1980; Snodgrass & Ulrich, 1990), скорость которой была получена по прямым допплеровским наблюдениям и быстрее внешней 5% области конвективной зоны (Beck & Schou, 2000), скорость которой была определена с помощью гелиосейсмологип. Еще более загадочным является тот факт, что супргра-нуляционная сетка вращается быстрее, чем магнитное поле (Snodgrass & Ulrich, 1990), причем скорость вращения зависит от размеров ячеек. Крупномасштабные структуры вращаются быстрее, чем мелкомасштабные (Duvall, 1980; Beck & Schou, 2000). Этот феномен был недавно подтвержден для поля скоростей, полученного по наблюдениям/-моды с помощыо локальной гелиосейсмологип (Gizon et al., 2003). В свете этих необычных свойств, а так же при отсутствии убедительных доказательств конвективного происхождения суп ер грануляции, для ее объяснения широко привлекаются альтернативные теории: модулирование конвективных движений гравитационными волнами (Lindzen & Tung, 1976), изменение с глубиной фактора заполнения магнитного поля (Foukal, 1977), взаимодействие конвекции и r-мод (Wolff, 1995), суперпозиция бегущих волн неизвестного происхождения (Gizon et al., 2003). В данной ситуации численное моделирование солнечной конвекции приобретает первоочередное значение.
Существует два дополняющих друг друга подхода к описанию солнечной конвекции. Первый использует упрощенную физику, чтобы исследовать свойства конвекции в глубоких слоях с учетом кориолисовых сил (Brummel ct al., 1996; Weiss et al., 1996; Elliott et ah, 1998). Как правило, эти расчеты невозможно продолжить до самой поверхности из-за введенных упрощающих предположений (Miesch, 2000). Другой подход, использующий реалистичные трехмер- ные расчеты конвекции (Атрощенко, 1993; Atroshchenko & Gadun, 1994; Stein & Nordlund, 1998, 2000), с большой точностью учитывает уравнение состояния вещества, перенос излучения в линиях, магнитное поле, и дает хорошее согласие с наблюдениями. Однако, такие расчеты требуют больших затрат процессорного времени и охватывают область всего в 1-2 гранулы. Для уверенного моделирования супергрануляции вычислительная область должна сдержать десятки гранул, находиться вблизи поверхности и учитывать сжимаемость среды.
Из-за непрозрачности поверхностных слоев Солнца, невозможно непосредственно наблюдать подфотосферные движения вещества, и на настоящее время, практически единственным инструментом для исследования внутреннего строения Солнца посредством прямых наблюдений является гелиосеисмология. Существует два различных взаимно дополняющих друг друга метода восстановления внутреннего строения по наблюдательным гелиосейсмологиче-ским данным. Первый подход основан на исследовании резонансных свойств внутренних областей Солнца путем определения частот нормальных мод колебаний (глобальная гелиосеисмология). Глобальная гелиосеисмология добилась впечатляющих успехов в исследовании крупномасштабных статических и динамических структур внутри Солнца (глубина конвективной зоны, дифференциальное вращение), что позволило значительно продвинуться в понимании физики процессов, протекающих внутри Солнца и других звезд. Другой подход основан на локальном исследовании бегущих волн (локальная гелиосеисмология), что позволяет исследовать явления меньшего масштаба, которые глобальная гелиосеисмология не в состоянии разрешить. Сюда входит исследование структуры солнечных пятен, магнитного поля и потоков вещества в супергрануляционных ячейках. Существует несколько методов исследования взаимодействия акустических волн с локальными неодиородиостями модели. Один из них — подход, использующий технику время-расстояние (time-distance), предложенную Дювалем с соавторами (Duvall et al., 1997). Другими альтернативными подходами являются анализ кольцевых диаграмм (ring-diagram analysis) и гелиосейсмическая голография (helioseismic holography). Ключевой концепцией подхода время-расстояние является понятие времени распространения волны между фиксированными точками на поверхности Солнца, которое можно получить путем вычисления кросс-корреляции колебаний на поверхности. Распространение звуковых колебаний при этом подходе вычисляется в рамках геометрической акустики. Однако, приближение геометрической акустики не применимо вблизи поверхности, где давление и плотность меняются очень быстро и должны приниматься во внимание волновые эффекты (Bogdan, 1997). В частности, время распространения звуковой волны чувствительно к свойствам среды не только вдоль пути распространения волны, но и к условиям в окрестности траектории. Более того, ядра чувствительности, вычисленные в борновском приближении (первое волновое приближение), имеют нулевое значение вдоль траектории (Kosovichev, et al., 2000). Таким образом, переход от геометрической акустики к волновой теории не сводится к простому уширению ядер чувствительности. Важным следствием является то, что чувствительность вдоль траектории распространения уже не пропорциональна обратной локальной скорости звука и может сильно от нее отличаться, особенно вблизи поверхности (Stark & Nikolaev, 1993).
В этой ситуации исключительно большую роль приобретает прямое численное моделирование взаимодействия распространяющихся звуковых волн с неоднородностями среды. Одномерные тестовые расчеты были выполнены Ко-совичевым и Дювалем (Kosovichev & Duvall, 1997), двумерные расчеты были проделаны Иенсеном, Якобсеном и Кристенсеном-Далсгартом (Jensen et al., 1999). Необходимо было провести реалистичные трехмерные расчеты рассеяния звуковых волн на неоднородностях внутреннего строения. Наблюдениями SOHO было подтверждено, что вспышки на поверхности Солнца могут порождать бегущие сейсмические волны и гигантские "солнцетрясения", подобные наблюдавшемуся 6 июля 1996 г.
Глобальная гелиосейсмология позволяет исследовать ис только крупномасштабные статичные структуры внутри Солнца, но так же глобальную дина- мику внутренних областей — дифференциальное вращение. В спектре мощности акустических "пятиминутных" мод колебаний наблюдается тонкая структура пиков, связанная с вращением Солнца (Claverie et at., 1981; Duval & Harvey, 1984; Brown & Morrow, 1987). Измеренные разности частот между пиками в тонкой структуре (так называемое вращательное расщепление) могут быть использованы для нахождения закона внутреннего вращения, который представляет большой интерес для изучения внутреннего строения и динамики Солнца, теорий звездной конвекции и механизма динамо, а также для определения гравитационного квадрупольного момента Солнца (Gough, 1985; Toomre, 1984; Макаров и др., 1987).
Смещения частот собственных колебаний вызваны, в основном, простым переносом волновой картины колебаний относительно наблюдателя в результате вращения. Величина этого смещения пропорциональна угловой скорости вещества внутри Солнца, усредненной с некоторым весовым множителем по области распространения колебаний. Размеры этой области по радиусу и широте различны для разных мод акустических колебаний. Поэтому, зная вращательное расщепление частот для некоторого набора мод, можно попытаться определить пространственную структуру внутреннего вращения. С математической точки зрения эта задача заключается в восстановлении зависимости угловой скорости от радиуса и широты по известным интегральным средним, то есть представляет собой обратную некорректную задачу для интегрального уравнения (Тихонов и Арсенин, 1979). В настоящее время разработаны методы решения этой задачи для определения радиальной зависимости угловой скорости в плоскости экватора (Duvall et al., 1984; Gough, 1984), а также широтного дифференциального вращения, заданного в параметрическом виде: или аналогичном представлении через полиномы Лежандра (Durney et al., 1987; Косовичев, 1988). Как известно (Вандакуров, 1967; Sh'ibahashi, 1979), высокочастотные акустические моды собственных колебаний могут быть рассмотрены в приближении ВКБ, поскольку длина волны этих колебаний вдоль радиуса ма- ла по сравнению с характерными масштабами изменения плотности и температуры. При этом для функций, описывающих зависимость угловой скорости от радиуса, получаются интегральные уравнения Абеля, решения которых находятся в явном виде по известной формуле абелевой инверсии (Gough, 1984). Однако, параметрическое представление П(л, 9) может не достаточно хорошо соответствовать реальному закону солнечного вращения. Например, угловая скорость в конвективной зоне может быть зависящей, в основном, от расстояния от оси вращения г sin 0 (Durney, 1976).
Ли и Шибахаши (Lee & Shibahashi, 1986) разработали метод определения 0(г, 0) в общем виде. Они предложили находить зависимость от угла, аппроксимируя ее кусочно-постоянной функцией и решая получающуюся при этом систему линейных уравнений. Однако, данная система может оказаться плохо обусловленной и потребовать задания дополнительной информации об угловой зависимости П, как это часто бывает при подобном подходе к численному решению интегральных уравнений I рода (Тихонов и Арсении, 1979). Мы рассмотрели моды колебаний с большими значениями угловой степени /, т. е. колебания с короткими длинами волн не только в радиальном, но и в поперечном направлении. Данное приближение представляет интерес по двум причинам: во-первых, акустические моды с /> 10 захвачены при г > 0.35 R , т. е. позволяют исследовать значительную область внутри Солнца; во-вторых, данные наблюдений, во время выполнения работы, были получены, в основном, для таких мод.
Кориолисовы силы, вызываемые вращением, должны оказывать существенное влияние на другие крупномасштабные движения вещества внутри Солнца, в частности на супергрануляциго. Непосредственные наблюдения супергрануляционной сетки затруднены, так как на относительно медленные движения вещества, связанные с суперфануляцией, накладывается быстрая хаотическая фануляционная картина. Суперфануляционная сетка становится видна только после соответствующей обработки поля горизонтальных скоростей отдельных фанул (local correlation tracking). В этой ситуации ключевую роль начинают играть методы восстановления и обработки изображений, позволяющие исследовать тонкую структуру солнечных образований: протуберанцев, пятен, волокон, грануляционной сетки. В данной работе предлагается метод определения аппаратной функции телескопа для последующего восстановления изображений солнечных протуберанцев. Наблюдения протуберанцев исключительно важны, так как они располагаются вдоль нейтральных линий магнитного поля, выявляя крупномасштабную картину, которая, в свою очередь, определяется глубинными конвективными структурами, которые исследуются с помощью гелиосейсмологии.
Влияние атмосферы и неидеалыюсти телескопа приводит к размыванию изображения и, как следствие, ухудшению разрешающей способности. В настоящее время используется несколько принципиально различных подходов к решению проблемы восстановления солнечных изображений.
Метод спекл-интерферометрии широко используется при исследовании мелкомасштабных солнечных структур. Метод позволяет восстанавливать изображение, вычисляя фазу истинного изображения из кросс-коррелянионного спектра, полученного из Фурье-образов последовательных кадров. Техника спекл-интерферометрии была предложена в работе (Labeyrie, 1970) и получила дальнейшее развитие в работах (Knox & Thompson, 1974; Weigelt & Wirnitzer 1983). Детальное описание реализации спекл-алгоритма, адаптированного специально для восстановления солнечных изображений дано в (von der Luhe, 1993). Для применения этого метода требуется серия («100) снимков одного и того же объекта с малой экспозицией. Метод позволяет свести к минимуму влияние атмосферы.
Метод фазового разнесения, являющийся разновидностью метода фокального объема. Для восстановления искажений волнового фронта используется информация о распределении интенсивности в трехмерной окрестности фокальной плоскости. Принципы метода фазового разнесения подробно описаны в работах (Gonsalves, 1982; Paxman et al., 1992; Lofdahl & Scharmer, 1994; Paxman et al., 1996). Использование метода применительно к наблюдению сол- печных пятен продемонстрировано в работе (Tritschler et al., 1997). На практике, как правило, получают два изображения: одно в фокальной плоскости, другое — слегка не в фокусе. Степень дефокусировки должна быть точно известна.
Метод интегральных уравнений (Тихонов и Арсенин, 1979) использует в качестве дополнительной информации значение аппаратной функции (АФ) К(ху), после чего решается интегральное уравнение типа свертки относительно истинного изображения. Наблюдаемое сглаженное изображение входит в правую часть уравнения.
Получающееся интегральное уравнение Фредгольма 1-го рода представляет собой типичную некорректно поставленную задачу, для решения которой необходимо применять те или иные регуляризируюшие алгоритмы. Понятие корректности ввел Ж. Адамар (Hadamard, 1902, 1932). Задача называется корректной по Адамару если 1) решение задачи существует; 2) решение единственно; 3) решение устойчиво. Последний пункт означает, что малые изменения исходных данных приводят к малым изменениям решения.
Для решения задачи восстановления изображения методом интегральных уравнений требуется одно изображение и знание АФ. На виде АФ сказывается не только неидеальность прибора, но так же и влияние атмосферы: неоднородности оптической плотности (атмосферные линзы) и рассеяние света. В идеале, АФ должна восстанавливаться из того же снимка, так как дрожание атмосферы будет приводить к тому, что АФ будет меняться от снимка к снимку. Это также важно при обработке архивных фотографий, когда невозможно провести дополнительные наблюдения для определения АФ. Как правило, достаточно трудно отделить вклад в АФ от самого прибора и атмосферы, но для восстановления изображения это не требуется, так как нас интересует суммарный эффект. В дальнейшем, под термином АФ мы будем понимать функцию, получающуюся при суммарном учете вклада атмосферы и неидеалыюсти телескопа.
Идея использования края Солнца для восстановления АФ возникла давно, например (Степанов, 1957). В предположении, что аппаратная функция имеет гауссову форму, в этой работе был предложен метод оценки параметра, опреде- ляющего полуширину АФ. Однако, на практике, вид АФ может сильно отличаться от гауссова (АФ может иметь уплощенную вершину или быть сильно скошенной). Использование гауссовой АФ вместо истинной приведет к большим ошибкам в восстановленном изображении.
При исследовании уравнения состояния и переноса излучения на Солнце важную роль играет поведение пылевой компоненты. Динамика пылевой компоненты начинает играть существенную роль при изучении звезд, находящихся на ранних этапах эволюции. Процесс седиментации (оседание пыли) также играет важную роль при исследовании формирования протопланетиых дисков. С другой стороны, оседание нерастворимой суспензии в вязкой среде играет важную роль в различных прикладных разделах физики, геологии, экологии. В промышленности седиментация часто используется для разделения суспензий на фракции, состоящие из частиц с различными радиусами. Таким образом, изучение седиментации полидисперсных суспензий представляет большой интерес. Спектр размеров частиц является важной характеристикой полидисперсной суспензии. Не смотря на то, что наиболее естественным способом исследования спектра размеров является построение гистограммы распределения частиц по некоторой выборке, позволяя оценить такие характеристики распределения как скошенность, мультимодальность, положения пиков и т.д., это редко делается. Как показано (Krumbein, 1934), такая гистограмма не является достоверной оценкой для распределения частиц, так как сильно зависит от размеров отобранных частиц. В той же статье предлагается использовать метод графического дифференцирования седиментационной кривой (массовой доли осевшего вещества суспензии в зависимости от времени) для определения функции распределения частиц суспензии по радиусам. Этот метод с небольшими вариациями (Шелудко, 1960) используется и в настоящее время. Обычная техника заключается в построении касательных к седиментационной кривой и измерения длин отрезков, ограниченных точками пересечения касательных и оси OY. Длины отрезков представляют собой массовую долю частиц суспензии в соответствующем диапазоне радиусов частиц. Такая методика позволяет получить только достаточно грубую гистограмму. Точки, в которых проводятся касательные к седиментационной кривой, должны быть достаточно далеко разнесены, так как большая погрешность в определении угла наклона касательных, проведенных в близких точках, даст большой флуктуационный вклад в высоту столбцов гистограммы. Соответственно, ширина отдельных столбцов гистограммы также должна быть достаточно большой. Необходим метод, свободный от вышеперечисленных недостатков. Большинство классических непараметрических методов оценивания плотности распределения так же не приводят к желаемому результату (Izenman, 1991). Для оценки плотности распределения частиц по радиусам больше подходят методы, основанные на сплайн-интерполяции (Wahba, 1975 a, b; Morandi & Costantini, 1989; Beatson & Wolko-wicz, 1989) или сплайн-аппроксимации (Вапник, 1984). Предлагаемый метод не опирается на аппроксимацию функции распределения сплайнами и сводится к задаче решения интегрального уравнения Вольтерры 1-го рода. В некотором смысле, эта задача занимает промежуточное положение между задачей решения уравнения Вольтеры П-го рода, которая всегда корректна и задачей решения уравнения Фредгольма 1-го рода, которая является некорректно поставленной в любых "разумных" функциональных пространствах. Задача решения уравнения Вольтерры 1-го рода может быть либо корректной, либо некорректной в зависимости от того, в каких пространствах она рассматривается (Апарцин, 1973, 1976, 1979, 1981; Апарцин и Бакушимский, 1972). Однако, даже в случае если задача корректна, некоторые численные алгоритмы могут порождать неустойчивости. Можно показать (Апарцин, 1973; Верлань и Сизиков, 1986), что методы, основанные на квадратурных формулах Сим пеона, Ныото-на-Котеса и других формулах высокого порядка, порождают расходящиеся алгоритмы. Формулы трапеций и прямоугольников приводят к устойчивым алгоритмам только в случае специального выбора шага интегрирования, зависящего от ошибки правой части и ядра. Таким образом, эти алгоритмы можно трактовать как регуляризирующие, в которых роль параметра регуляризации играет шаг интегрирования. Это накладывает определенные ограничения на использо- вание этих методов, так как правая часть, как правило, известна из эксперимента в виде таблицы значений на нерегулярной сетке. Таким образом, для получения значений правой части в узлах равномерной сетки с оптимальным шагом, нам необходимо решить нетривиальную (в общем случае) задачу интерполяции зашумленных экспериментальных данных.
Уравнение Вольтерры 1-го рода можно формально рассматривать как специальный случай уравнения Фредгольма 1-го рода с ядром, содержащим ступенчатую функцию. Следовательно, регуляризирующие алгоритмы, разработанные для решения уравнений Фредгольма 1-го рода, с успехом могут быть применены и к решению интегральных уравнений Вольтерры 1-го рода. Здесь всюду в дальнейшем применялся метод регуляризации А.Н. Тихонова.
Рассмотренные методы восстановления функции распределения частиц по радиусам, применимы в случае, когда суспензия оседает в неподвижной среде. Очень часто приходится рассматривать поведение суспензии в движущейся среде. В некоторых случаях движение среды имеет выделенное направление, как при течении жидкости по трубам, в других, как в случае конвекции, — нет. Конвективные движения среды возникают из-за наличия температурных градиентов. Если мы не предпринимаем специальных мер по стабилизации температуры, сплошная среда почти всегда приходит в конвективное движение. В этой ситуации исключительно важно знать какое влияние окажет конвекция на процесс седиментации, и как исказятся результаты восстановления функции распределения частиц по радиусам в присутствии конвекции. Херман (Herrmann et al., 1999) провел численное моделирование седиментации монодисперсной и бидисперсных суспензий в статичной вязкой среде, однако, эти результаты не могут быть прямо перенесены на седиментацию полидиспресной суспензии в конвективно-неустойчивой среде. Для этого требуется дополнительное моделирование. Большой интерес представляет сравнение картины конвективных движений в субфотосферных слоях Солнца и в лабораторных условиях. Выявляются принципиальные качественные различия, связанные с большим градиентом плотности вещества на Солнце и большим значением ускорения свобод- ного падения. На Солнце характерной особенностью конвекции является наличие пространственно узко локализованных высокоскоростных нисходящих потоков холодного вещества, которое у поверхности быстро остывает как за счет адиабатического расширении при подъеме, так и за счет излучения. Скорости поднимающегося на большой площади горячего вещества существенно ниже. В лабораторных условиях картина противоположная. Горячее поднимающееся вещество собирается в высокоскоростной тонкий "жгут", в то время, как холодное вещество медленно опускается на большой площади.
Очень часто при реализации алгоритмов, предлагаемых в данной работе, нам приходилось численно вычислять производную от зашумленных экспериментальных данных. Эта задача имеет и самостоятельное значение. Понятие "скорости процесса", например, исключительно важно в биологии, где оно играет фундаментальную роль. Помимо очевидного использования понятия скорости при исследовании роста организмов, существует масса других областей применения. Жизнедеятельность любого отдельного организма связана с дыханием, питанием, ростом, (у растений добавляется фотосинтез, при исследовании целых экосистем необходимо принимать во внимание миграцию особей, изменение численности популяции). Все эти процессы могут быть описаны в терминах потоков кислорода, питательных веществ, скоростей роста, миграции. Так или иначе, но для количественного описания этих процессов необходимо уметь измерять (или вычислять) скорости. Как правило, в эксперименте не измеряется непосредственно скорость процесса, ее приходится вычислять путем численного дифференцирования.
На фоне важности понятия «скорости» для исследования биологических объектов, не понятен разнобой, присутствующий в биологической литературе в использовании этого термина. Помимо мгновенной скорости и удельной скорости v = V v Л& dt' "' f d/' =7^7 !/<'>*. часто используются средние скорости на некотором интервале времени. Операция временного усреднения () некоторой функцииУ(/) по временному интервалу [t\, /2] определяется следующим образом:
Используя эти определения нетрудно получить выражения для средних скоростей в виде
If,"' 4 4 ^in^l
1 Uw А. /a-'l /(',)
Необходимо подчеркнуть, что выражения для средних скоростей, справедливы для произвольной зависимости/от времени. Очень часто, в биологической литературе авторы опускают слово "средняя" в понятии "средняя скорость", что может приводить к недоразумениям. Только в случае, когда/зависит от t линейно, средняя скорость совпадаете мгновенной скоростью. Аналогично, средняя удельная скорость совпадает с мгновенной удельной скоростью только для случая экспоненциального поведения функции J(t). Второе уравнение было впервые введено в биологии в 1927 году одновременно американским и русским учеными Броди и Шмальгаузеном (Brady, 1945; Шмальгаузен, 1984а,б). Они рассматривали его как мгновенную удельную скорость роста, неявно предполагая экспоненциальный закон роста.
Наряду с неточностями в использовании понятия "скорость", присутствуют явные ошибки. Например, множество исследователей (DeBoer et al., 1978; Fujita & Goldman, 1985; Gordon et al., 1981; Hwang et al., 1987; Lobban et al., 1985; Liming, 1990), следуя друг другу, умножают (Vfci) на 100% и интерпретируют получающуюся величину как скорость увеличения массы, выраженной в процентах, не замечая, что Vre\ — размерная величина. Некоторые исследователи (Винберг, 1968; Карманова, 1976; Neish et al., 1977), для удельной средней скорости, используют выражение w"'y /,-/, /0,) или его модификацию w- /0,)-Л',) ':-' ][/0,) + /02)]* it Эти же выражения используются авторами (Хайлов и др., 1992; Lebedev et al., 1989) для исследования потоков веществ в экосистемах. Эти выражения не только не дают правильное значение средней удельной средней скорости, они вообще не являются средними какой-либо величины на произвольном интервале [/i, t2]. Некоторые исследователи (Haglund & Pedersen, 1992, 1993; Lignel et al., 1987; Ohno &. Mairth, 1982; Soe-Mtun et al., 1986) используют выражение 100%, V.- = (filA) '" I/0,) J ошибочно трактуя его как скорость роста. На самом деле, величины Vre{ и К% связаны следующим соотношением: х г,,/ 400% J Такое разнообразие ошибочного использования понятия "скорость" в биологической литературе и побудило нас написать небольшое введение (Парчевский и Парчевский, 1998). Средняя скорость является хорошим интегральным показателем, который хорошо описывает процесс в целом. Однако, если мы начинаем дробить интервалы усреднения на более мелкие для того, чтобы проследить динамику процесса, ошибки средней скорости начинают увеличиваться. Ошибки могут стать в несколько раз больше самой величины и мы, вместо гладкой функции, получим пилообразную функцию большой амплитуды. Таким образом, мы видим, что средняя скорость не подходит для изучения динамики процессов с достаточно хорошим временным разрешением. Нам необходимо использовать другие подходы для восстановления производной из экспериментальных данных. В данной работе приводится сравнение подходов, основанных на сплайн-аппроксимации и решении интегрального уравнения на производную методом тихоновской регуляризации. С задачей восстановления производной тесно связана другая не менее важная задача — восстановление функции плотности вероятности (ФПВ) по выборке конечного объема. Стандартный метод разбиения на классовые интервалы с последующим подсчетом количества точек в каждом интервале приводит к грубой гистограмме, не дающей представления о тонкой структуре ФПВ. Чтобы получить непрерывную ФПВ необходимо аппроксимировать полученную гистограмму каким-либо выборочным распределением и подобрать его коэффициенты с использованием критерия у}. Как правило, малое количество свободных параметров в выборочных распределениях не дает возможности детально описать тонкую структуру ФПВ. В случае сложных многомодовых распределений вообще невозможно подобрать подходящую выборочную функцию плотности. Большинство классических непарамстрических методов оценивания плотности распределения не приводят к желаемому результату (Izenman, 1991), Один из подходов заключается в использовании аппроксимирующих (Вапник, 1984) или интерполирующих (Wahba, 1975 a, b; Morandi & Costantini, 1989; Beatson & Wolkowicz, 1989) кубических сплайнов. Другой подход основан на решении интегрального уравнения, связывающего кумулятивную функцию распределения и функцию плотности, и, по-сути, сводится к вычислению производной от эмпирической кумулятивной функции распределения. Метод позволяет получать непрерывную гладкую ФПВ и свободен от недостатков методов, использующих сплайн-аппроксимацию. Цели и задачи исследования. Используя упрощенную физику, провести моделирование сжимаемой солнечной конвекции в двух измерениях в области, охватывающей по горизонтали несколько десятков гранул и простирающейся до глубины, где высокотурбулентные высокоскоростные приповерхностные движения вещества будут сильно подавлены. Основываясь на результатах численного моделирования, исследовать эволюцию отдельных гранул, их взаимодействие и гори- зонтальное движение, выявить наличие крупномасштабных структур, соответствующих на Солнце мезо- и супергрануляции. Используя реалистичную физику, промоделировать в трех измерениях рассеяние звуковых волн на солнечном пятне. Исследовать задержку и изменение амплитуды волнового фронта при прохождении сейсмической волны через солнечное пятно. Построить численную модель "солнцетрясения" — распространения сейсмической волны от солнечной вспышки, наблюдавшейся 6 июля 1996 г. Исследовать численную модель на наличие колебаний и волн, возбуждаемых конвекцией и найти их спектр. Провести дальнейшее развитие асимптотического метода для определения внутреннего дифференциального вращения и разработать непараметрический метод для определения зависимости угловой скорости дифференциального вращения Солнца от широты и глубины. Разработать непараметрический метод, позволяющий восстанавливать АФ из наблюдений солнечного лимба, причем АФ должна восстанавливаться из того же снимка, который в дальнейшем будет подвергнут обработке с целью увеличения разрешения. Разработать непараметрический метод, позволяющего восстановить непрерывную гладкую функцию распределения частиц суспензии по радиусам, используя экспериментальную седиментационную кривую. Провести численное моделирование седиментации полидисперсной суспензии в конвективно-неустойчивой среде и выявить характерные закономерности поведения такой системы. Результаты, выносимые на защиту. 1, Проведенное моделирование двумерной сжимаемой солнечной конвекции доказывает существование спонтанно возникающих крупномасштабных структур, по размеру и времени жизни соответствующих супергрануляциопным ячейкам на Солнце. Показано, что турбулентная конвекция вблизи дна вычислительной области возбуждает гравитационные волны с дискретным спектром. Построена численная модель "солнцетрясения" — распространения сейсмической волны от солнечной вспышки, наблюдавшейся спутником SOHO 6 июля 1996 г. Расчет рассеяния акустических волн на солнечном пятне позволил определить задержку волнового фронта рассеянной волны. Показано, что амплитуда звуковых волн в пятне увеличивается, поверхностная скорость также увеличивается по мере удаления фронта волны от источника. Предлагается метод восстановления угловой скорости внутреннего дифференциального вращения Солнца по наблюдаемому вращательному расщеплению 5-мин. колебаний. Показано, что в асимптотическом приближении задача может быть сведена к решению двумерного интегрального уравнения Абеля. На угловую скорость вращения априори не накладывается никаких функциональных ограничений. Получено аналитическое решение задачи в квадратурах. Предлагается непараметрический метод определения аппаратной функции солнечного телескопа для восстановления изображений протуберанцев. Аппаратная функция определяется из На наблюдений спокойного края солнечного диска как решение двумерного интегрального уравнения типа свертки. Показано, что сечение двумерной аппаратной функции плоскостью X0Z равно производной по х от профиля яркости, усредненного вдоль лимба. Предлагается новый непараметрический метод восстановления функции распределения частиц полидисперсной суспензии по радиусам. Подход основан на решении интегрального уравнения, описывающего процесс седиментации (оседания) суспензии в вязкой среде, методом тихоновской регуляризации. Предлагаемый метод позволяет получать гладкую функцию распределения сложной многомодовой формы и может использоваться для анализа смеси двух и более суспензий с различными функциями распределения частиц по размерам. Проведенное численное моделирование седиментации полидисперсной суспензии в прямоугольной области в конвективно-неустойчивой среде показало, что конвекция действует как размерный фильтр, сепарируя частицы суспен- зии на взвешенную и осевшую фракции в зависимости от их радиусов. Радиусы частиц во фракциях зависят только от разности температур верхней и нижней граней. Научная новизна полученных результатов. Численное моделирование двумерной сжимаемой солнечной конвекции показало, что супергрануляционные структуры спонтанно возникают при взаимодействии соседних гранул. Впервые построена численная модель события 6 июля 1996 г. — распространения сейсмической волны от солнечной вспышки, наблюдавшейся спутником SOHO. Впервые проведено реалистичное 3D моделирование рассеяния акустических волн солнечным пятном. Подход к обработке различных наблюдательных и экспериментальных данных с единой позиции решения обратных некорректных задач позволяет разработать новые непараметрические методы для 1) восстановления угловой скорости внутреннего дифференциального вращения Солнца, 2) определения аппаратной функции солнечного телескопа. 3) восстановления функции распределения частиц полидисперсной суспензии по радиусам. Численное моделирование седиментации полидисперсной суспензии в прямоугольной области в конвективно-неустойчивой среде показало, что конвекция действует как размерный фильтр, сепарируя частицы суспензии на взвешенную и осевшую фракции в зависимости от их радиусов. Личный вклад соискателя. Диссертационная работа является самостоятельным научным исследованием. Лично автором разработан план исследований. Используемые алгоритмы распараллелены и реализованы в виде пакетов программ на C++. Для расчета солнечной конвекции и распространения акустических волн были отдельно разработаны распараллеленные версии программ, оптимизированные для работы на векторных многопроцессорных системах. Программа восстановления изображений солнечных протуберанцев реализована в среде Matlab. При вычислении коэффициентов численной схемы при расчете распространения акустических волн использовалась система аналитических вычислений Mathematica. Расчеты проводились на суперкомпьютерах Отдела суперкомпьютерных вычислений исследовательского центра НАСЛ, Центра турбулентных исследований Стэнфордского университета и суперкомпьютерах Ханзеновской лаборатории экспериментальной физики Стэнфордского университета (Калифорния, США). Лично автором проведена интерпретация полученных результатов. Апробация работы. Материалы диссертации докладывались автором на семинарах Отдела физики Солнца Крымской астрофизической обсерватории (НИИ КрАО) 1990-2003 гг., заседании президиума Академии наук Украины (Киев, 1999 г.), конференциях по физике Солнца (НИИ КрАО, 1998, 1999, 2001, 2002, 2003 гг.), первой международной конференции по математической экологии (Мадрид, Испания, 3-8 сентября 1998 г.), рабочей группе по математическим проблемам потоков суспензий (Штутгарт, Германия, 9-Ю октября 1999 г.), конференции по сепараторным системам жидкость-твердое тело (Гавайи, США, 18-23 апреля 1999 г.), конференции по эволюционному моделированию систем частиц (Гавайи, США, 23-28 января 2000 г.), рабочей группе по локальной гелиосейсмологии (Стэнфордскин университет, Калифорния, США, 6-7 апреля 2001 г.), конференции по локальной и глобальной гелиосейсмологии SOH012/GONG+ (обсерватория Биг Бэр, США, 27 октября - 1 ноября 2002 г), объединенном семинаре Института астрономии Кембриджского университета (ГОА, Кембридж, Великобритания), семинарах Исследовательского центра НА-СА (NASA Ames Research Center, Калифорния, США), Ханзеновской лаборатории экспериментальной физики Стэнфордского университета (HEPL, Стэн-фордский университет, Калифорния, США), Центра турбулентных исследований Стэнфордского университета (CTR, Стэнфордский университет, Калифорния, США), Геологоразведочной службы США (US Geological Survey, Калифорния, США). Публикации. Результаты диссертационной работы полностью опубликованы в 15 научных работах (6 без соавторов), в числе которых 8 статей в зарубежных изданиях. Содержание диссертации. Диссертация состоит из введения, трех глав и заключения. Во введении обосновывается актуальность темы, дается литературный обзор по выделенным проблемам и, исходя из этого, проводится постановка задач. Первая глава посвящена исследованию крупномасштабных свойств солнечной конвекции и рассеянию акустических волн на неоднородностях внутреннего строения Солнца путем численного моделирования. Во второй и третьей главах разрабатываются новые методы обработки различных наблюдательных и экспериментальных данных с единой позиции решения обратных некорректных задач. Вторая глава посвящена разработке новых методов обработки солнечных наблюдений: восстановлению внутреннего дифференциального вращения из гелиосейсмологических данных, определению аппаратной функции телескопа по наблюдениям солнечного лимба и восстановлению Нц изображений протуберанцев. В третьей главе разрабатываются методы решения различных прикладных задач промышленности и экологии, имеющих большое практическое значение: определение размеров частиц суспензий, разработка методов тонкого фракционирования суспензий с использованием конвекции, определение скорости роста организмов. Моделирование солнечной сжимаемой конвекции проводилось в прямоугольной области с отношением сторон 6:1 и сеткой, содержащей 600x100 узлов. Геометрический размер вычислительной области составлял L = 47.95 Мм по горизонтали и Н 7.99 Мм по вертикали. Отношение плотностей среды на нижней и верхней гранях равнялось рь/р, = 1800. Для интегрирования системы уравнений (6) в безразмерных переменных использовался явный метод Мак Кормака (12, 13). Шаг по времени Af = 1.04 с выбирался в соответствии с ограничением Куранта (16). Общая продолжительность моделирования составляла 9.317 суток солнечного времени (775000 итераций). После окончания переходного процесса на 265000-й итерации начали выводиться результаты счета. Промежуточные результаты сохранялись с интервалом в 100 итераций (A/out = 104 с) в виде файлов, содержащих значения температуры, плотности и скоростей, а также в виде фильма, кадрами которого являлись карты распределения температуры. Фильм был представлен в Отделе суперкомпьютерных вычислений Центра турбулентных исследований в Исследовательском центре НАСЛ. Фильм в формате MPEG, можно найти на сайте На рис. 1 показан один из кадров фильма. Яркостью помечена разность температур среды нашей модели и модели с линейным ростом температуры. Темные участки обозначают области с температурой ниже, чем в модели с постоянным градиентом. Светлые — с более высокой температурой. Стрелками на рис. 1а помечен полный вектор скорости среды. Сохранено соотношение 6:1 сторон вычислительной области. На рис. 1б,в показаны, соответственно, карты вертикальной и горизонтальной составляющих скорости. Характерной особенностью солнечной конвекции является быстрое остывание вещества в приповерхностном слое и образование, так называемых, да-ундрафтсов (в иностранной литературе downdrafts) — компактных, сравнительно холодных образований, движущихся вниз с большой скоростью 10.433 ± 0.046 км/с и проникающих глубоко в область с более спокойными движениями. Они образуют тонкие стенки гранул. В центре гранулы происходит подъем горячего вещества со значительно меньшими скоростями 5.270 ± 0.023 км/с. Характерное время жизни гранулы составляет 10 -f- 15 мин. Чтобы проследить закономерности в движении гранул, необходимо построить так называемую диаграмму время-расстояние (time-distance diagram). Для этого выберем горизонтальный срез вычислительной области на интересующей нас глубине и проследим его эволюцию во времени. На рис. 2 показаны диаграммы время-расстояние температуры и скорости для слоя на глубине 0.55943 Мм, считая от верхней границы вычислительной области. По оси ОХ отложена горизонтальная координата в мегамстрах (1000 км), по вертикальной оси — время в часах. Нетрудно отождествить отдельные характерные детали в виде «елочек» на всех рисунках. На рис. 2а показано поле температуры. Темные области соответствуют холодному веществу, движущемуся вниз. Об этом говорит темная окраска соответствующих деталей на рис. 2в, представляющем карту вертикальной составляющей скорости. Характерные «елочки» на рис. 2а образуются благодаря горизонтальному движению даунд-рафтсов (стенок гранул), состоящих из холодного вещества. «Ветки», имеющие положительный наклон, формируются даундрафтсами, движущимися слева на право к центральному «стволу», отрицательный — даундрафтсами, движущимися справа налево. Это иллюстрируется рис. 2г на котором показана карта горизонтальной составляющей скорости. Левые «ветки» имеют светлую окраску, что говорит о том, что вещество там движется в положительном направлении оси ОХ (т.е. вправо). Правые «ветки», соответственно, имеют темную окраску. Там вещество движется влево. Горизонтальные движения даундрафтсов возникают из-за наложения на сравнительно мелкомасштабную высокотурбулентную картину грануляции более медленных крупномасштабных движений. При этом, поднимающееся горячее вещество у поверхности начинает растекаться в стороны, «смывая» как целое даундрафтсы вправо и влево от оси подъема. Вновь появляющиеся в этом месте даундрафтсы «смываются» как целое аналогичным образом, причем с той же горизонтальной скоростью, которая определяется горизонтальной скоростью растекания поднимающегося горячего вещества. Так формируются почти параллельные «ветки» у «елочек». Характерная долгожн-вущая «елочка» расположена на рис. 2 в окрестности 34 Мм. На рис. 26 показана диаграмма время-расстояние для модуля скорости среды. Здесь «елочки» имеют светлую окраску, в отличие от темных промежутков между ними. Это говорит о том, что скорости даундрафтсов, движущихся вниз, намного больше скоростей поднимающегося горячего вещества. Наблюдения показывают наличие на поверхности более крупномасштабных чем грануляция структур — супер- и мезогрануляции. Эти структуры не видны при непосредственных наблюдениях так как амплитуда соответствующих скоростей не высока и картина «замывается» высокоскоростными мелкомасштабными грануляционными движениями. Супергрануляция проявляется только при построении карт горизонтальной скорости гранул или диаграмм время-расстояние. Чтобы проследить возникают ли подобные крупномасштабные структуры при численном моделировании, мы построили диаграмму время-расстояние для всей длины вычислительной области L = 47.95 Мм, охватывающую интервал времени в 147.3 часа. Соответствующая иллюстрация приведена на рис. 3. На диаграмме время-расстояние, построенной для всей длины вычислительной области L = 47.95 Мм и охватывающей интервал времени в 147.3 часа, видны крупномасштабные структуры с характерным размером 15+20 Мм и характерным временем жизни 20 -т- 30 часов. Полученные результаты хорошо согласуются с супергрануляционной картиной, полученной из наблюдений SOHO/MDI (Shine et а!., 2000). Высокотурбулентная приповерхностная конвекция возбуждает различные типы волн: звуковые, гравитационные. Для поиска стоячих волн в нашей модели и определения для них дисперсионного соотношения было построено двумерное Фурье-преобразование диаграммы время-расстояние аналогичной рис. 4, но вычисленной для слоя па глубине h = 6.234 Мм, т.е. вблизи дна вычислительной области. Отчетливо просматриваются три характерных гребня подобные тем, что видны на спектре 5-мин. солнечных осцилляции. Работа опубликована в статьях (Parchevsky, 2001b, 2002b). Таким образом, в случае / » I, когда справедливо асимптотическое представление (34, 35) для присоединенных функций Лежандра, обратная задача восстановления дифференциального вращения Солнца по спектру наблюдаемых частот решается в квадратурах. Данное решение обладает существенными преимуществами по сравнению с другими решениями уравнения (36). В частности, в работе Ли и Шиба-хаши (Lee & Shibahashi, 1986) предлагается применить однократно инверсию Абеля по переменной у к уравнению (36), а оставшийся интеграл по в разбить на / + т интервалов и решать получающуюся систему линейных уравнений. Как уже отмечалось, этот метод неустойчив. Существенно проще и точнее выполнить численное интегрирование, если решение известно в квадратурах. При этом, однако, необходимо учитывать, что вследствие ошибок наблюдений для устойчивости инверсии Абеля требуются специальные регуляризируюшие алгоритмы (Войскобойников и др., 1984; Белоцерковский и Лифанов, 1985). Результаты. В данной работе предлагается метод восстановления угловой скорости внутреннего дифференциального вращения Солнца по наблюдаемому вращательному расщеплению спектра 5-минутных колебаний, причем априори не делается никаких предположений о функциональной зависимости Q(r, 0) от координат. В частности, метод применим даже в случае, если Q нельзя представить в виде произведения функций, одна из которых зависит только от г, а другая — только от в (то есть когда переменные не разделяются). Основная идея метода заключается в замене точных решений невозмущенной задачи их асимптотическим квазиклассическим приближением не только по переменной / , но и по переменной в. Это разложение справедливо в области /:» 1. В результате мы получаем для функции О. двумерное интегральное уравнение Абеля, правая часть которого зависит от двух автомодельных переменных: ,v =(/ + 1/2)/ст и и = (/ + \12)1т. Не составляет труда решить это уравнение в квадратуpax и получить явную зависимость Q(r, 0). При обработке наблюдений вращательное расщепление частот следует представлять в виде: о піт = щАх w)- Предложенный алгоритм был использован Секи (Sekii & Gough, 1993) для восстановления внутреннего дифференциального вращения Солнца. Результаты опубликованы в статье (Косовичев и Парчевскнй, 1988) и восстановление изображений солнечных протуберанцев. Предположим, что мы наблюдаем точечный объект. Идеальный телескоп дал бы нам точечное изображение. Если же мы используем реальный инструмент, мы получим не точку, а некоторую неравномерно освещенную область (см. рис. 10). Идеальная точка размывается. Функция, которая описывает распределение яркости, получающееся при наблюдении точечного источника, называется аппаратной функцией (ЛФ) прибора (в англоязычной литературе принят термин PSF - Point Spread Function). В случае, когда изображение содержит более одной точки, такому размыванию подвергается каждая точка исходного изображения, подобно сглаживанию с помощью бегущего среднего в двух измерениях. Такое инструментальное сглаживание приводит к ухудшению качества изображения и уменьшению разрешающей способности прибора. Нетрудно написать уравнение, связывающее исходное g( ,rf) и экспериментально полученное (сглаженное) Длу) изображения: Правая часть (40) представляет собой наблюдаемое искаженное и зашумленное изображение. Ядро К(ху) - есть АФ, которая может быть задана как аналитически, так и в виде некоторого двумерного массива экспериментальных точек, и также может содержать ошибки. Первый знак равенства можно рассматривать как определение интегрального оператора, переводящего истинное изображе-ниє в наблюдаемое. Выражение (40) может рассматриваться как интегральное уравнение на функцию (,?;), если нам известна, в том или ином виде, АФ, которая может быть получена экспериментально путем наблюдения точечного объекта. С другой стороны, если нам известно исходное изображение g(,T}), то (40) можно рассматривать как интегральное уравнение для определения ЛФ. Алгоритм восстановления изображений, основанный на решении интегрального уравнения (40), существенным образом использует информацию об ЛФ прибора. При реальных наблюдениях, как правило, АФ не известна заранее. В случае звездных телескопов не возникает проблем в получении аппаратной функции — достаточно получить изображение звезды. Из-за большой удаленности звезд это изображение с высокой точностью можно считать изображением точечного источника. По определению, наблюдаемое изображение и есть АФ телескопа. Ситуация существенным образом изменяется при наблюдениях Солнца. Как правило, на изображениях Солнца невозможно выделить точечные источники. Лабораторные источники света не являются точечными с необходимой точностью, а для наблюдений звезд солнечные телескопы, с их регистрирующей частью, не приспособлены. Для получения изображений протуберанцев и пятен используется узкополосный интерференционно-поляризационный На-фильтр, пропускающий менее 0.1% падающего на него света, что совершенно недостаточно для наблюдения даже самых ярких звезд. В данной работе предлагается непараметрический метод восстановления АФ солнечного телескопа по наблюдаемому изображению края солнечного диска. Метод основан на решении интегрального уравнения (40) относительно аппаратной функции К(х,,у) при известных g(trf) пДх у). Предположим, что на фотографии присутствует спокойный участок солнечного лимба, свободный от факелов и протуберанцев. Выделим этот участок и будем его использовать для получения АФ. При наблюдениях в крыле линии Ни, изображение края диска Солнца имеет очень резкую границу. Координаты и скорости частиц суспензии вычислялись на каждом временном шаге (12, ІЗ). В процессе численного эксперимента фиксировалась масса вещества, осевшего на нижней грани вычислительной области. Частица считалась осевшей, если на очередном временном шаге ее -координата оказывалась меньше или равна нулю. При этом ее масса добавлялась к массе осевшего вещества, а сама частица исключалась из дальнейших расчетов. Таким образом, численное моделирование дает нам не только седиментационную кривую P(t) (массовую долю осевшего вещества в зависимости от времени) но и информацию о радиусах частиц в осевшей и взвешенной фракциях суспензии. Для исследования распределения частиц по радиусам в каждой из фракций необходимо уметь восстанавливать функцию плотности вероятности по выборке конечного объема. Обычный метод восстановления функции плотности вероятности q{x) случайной величины х по выборке (х\, Хг ..., XN) конечного объема N состоит в следующем. Весь интервал изменения х разбивается на несколько подинтерва-лов, называемых классовыми интервалами. Далее, элементы выборки группируются в классы по признаку попадания в тот или иной классовый интервал и подсчитывается число элементов выборки щ (групповая частота) в j-ом подин-тервале. Результаты могут быть представлены графически в виде гистограммы. Величина классовых интервалов определяется исследователем. Если длина классовых интервалов мала, то в / будут доминировать случайные вариации, так как классовые интервалы будут содержать небольшое количество точек. При выборе больших классовых интервалов будет происходить сильное сглаживание функции плотности вероятности из-за усреднения. В реальных уеловиях выборок малого объема не всегда можно выбрать оптимальное разбиение. Кроме этого, данный подход не дает возможность получить гладкую функцию плотности вероятности. Для получения гладкой функции, обычно производятся следующие действия. Выбирается некоторая пробная функция плотности вероятности, и с ее помощью вычисляются относительные частоты nJN попадания в j-й классовый интервал. Теперь можно сравнить пробное распределение величины х с ее эмпирическим распределением. Используя критерий yj, можно либо принять, либо отвергнуть статистическую гипотезу о том, что пробное распределение совпадает с заданной достоверностью с эмпирическим распределением. Если статистическая гипотеза отвергается, мы должны выбрать другое пробное распределение и повторять процедуру до тех пор, пока статистическая гипотеза не будет принята. В случае сложных (многогорбых) функций распределения плотности вероятности, мы можем и не подобрать пробного распределения, которое с заданной достоверностью описывало бы наше эмпирическое распределение. Однако существует подход, позволяющий получить гладкую функцию плотности распределения вероятности из выборки конечного объема. Эта проблема тесно связана с задачей восстановления производной из зашумленных данных. Функция распределения F(x) непрерывной случайной величины х связана следующим образом с плотностью распределения вероятности q(x): Предположим, у нас есть упорядоченная по возрастанию выборка конечного объема (JCJ, Х2, ..., хм). Построим по данной выборке эмпирическую функцию распределения Fx) Функция Fe(x) есть ступенчатая неубывающая функция, являющаяся несмещенной состоятельной оценкой для функции распределения F(x). Чтобы найти функцию плотности вероятности q(x) необходимо подставить (65) в правую часть (64) и решить получившееся интегральное уравнение относительно функции q(x). Не трудно заметить, сравнивая (64) и (45), что задача восстановления плотности вероятности сводится к задаче нахождения производной от зашум-ленной экспериментальной функции. Обе задачи сводятся к однотипным интегральным уравнениям и являются некорректно поставленными. Интегральное уравнение (64) может быть решено различными методами. Мы использовали метод регуляризации Тихонова. Этот метод дает хорошие результаты и для выборок малого объема, когда правая часть (64) известна с большими ошибками. Таким образом, вычислительная процедура состоит из трех основных шагов: 1) численное моделирование конвекции до достижения системой квазистационарного состояния; 2) численное моделирование седиментации суспензии; 3) восстановление функции распределения частиц по радиусам в осевшей и взвешенной фракциях. В данном разделе исследовалась седиментация цементной пыли (р = 2200 кг/м3) в воздухе при нормальных условиях и наличии конвективных движений среды, вызванных градиентом температуры. Квадратная вычислительная область со стороной L - Юм содержала 50x50 равномерно расположенных узлов сетки. Состояние среды на верхней границе вычислительной области (опорное состояние) было выбрано следующим образом: Сама конвекция вызывается температурным градиентом между верхней и нижней гранями вычислительной области. Моделирование проводилось для двух случаев: AT = 20 К и AT = 50 К. Шаг по времени выбирался равным Д/ = 1.53 10"4 сек в соответствии с ограничением (16). В качестве граничных условий на верхней и нижней гранях использовалось условие глухой стенки: температура фиксирована, жидкость прилипает к стенке (скорости равны нулю). На боковых гранях ставились периодические граничные условия. В качестве начальны условий использовалось гидростатически равновесное состояние с постоянным градиентом температуры (температура линейно уменьшается от нижней грани до верхней). Для этого модельного случая мы уже получили аналитическое решение (18), где Г0, 7 1 играют роль температур нижней и верхней граней соответственно, Ду = L. Кривые представляют собой массовую долю осевшего вещества в зависимости от времени (время приводится в нормализованных единицах). Асимптоты (горизонтальные пунктирные линии) вычислялись следующим образом. Строилась кривая P{\lt) и экстраполировалась до пересечения с осью OY. Координата точки пересечения дает нам требуемую величину. Конвективные движения среды приводят к следующему эффекту. Мелкодисперсная фракция суспензии остается во взвешенном состоянии значительно дольше, чем в отсутствие конвекции. Некоторые частицы с достаточно малыми радиусами не оседают совсем. Массовая доля взвешенной фракции зависит от средних конвективных скоростей среды. Те, в свою очередь, зависят только от разности температур верхней и нижней граней. Чем больше эта разность, тем больше конвективные скорости и тем больше массовая доля взвешенной фракции. В данной ситуации конвекция действует как размерный фильтр, разделяя частицы на фракции в соответствии с их радиусами. Средний, а так же наиболее вероятный радиусы частиц во взвешенной и осевшей фракциях зависят от параметров конвекции. Эти параметры, в свою очередь, зависят только от разности температур верхней и нижней граней. Этот эффект дает возможность использовать конвекцию для тонкого фракционирования суспензий. После оседания крупных частиц они удаляются, после чего выключается подогрев, конвекция прекращается и взвешенная фракция оседает. Критический радиус частин, при котором происходит разделение на фракции, легко регулировать, изменяя разность температур между верхней и нижней гранями. Распределение частиц по размерам во взвешенной и осевшей фракциях представляет большой практический интерес. Так как во время моделирования фиксируются радиусы осевших частиц, то у нас есть выборки радиусов частиц взвешенной и осевшей фракций, и мы можем восстановить функции распределения частиц по размерам независимо для каждой из фракций. Результаты такого восстановления показаны на рис. 20а для АТ= 50 К и рис. 206 для АТ 20 К. Пунктиром показано исходное распределение частиц по радиусам (63). Более узкая кривая распределения, у которой максимум смещен в сторону малых радиусов (по отношению к исходному распределению) представляет собой распределение частиц во взвешенной фракции. Более широкая кривая распределения с максимумом, смещенным в сторону больших радиусов, дает нам распределение частиц в осевшей фракции. Методика восстановления функции распределения частиц по радиусам из седиментационной кривой, описанная в разделе 2.3 предполагает, что среда, в которой оседают частицы, неподвижна. Наличие конвективных движений среды искажает результаты, однако такие расчеты могут оказаться полезными. Результаты применения этой техники к седиментационным кривым на рис. 19 иллюстрируются рис. 21 и рис. 22. Пунктирной линией, как и на предыдущих рисунках, показано исходное логнормальное распределение (63). Сплошной линией на рис. 21 показана функция распределения частиц по радиусам, восстановленная по седиментационной кривой для АТ= 50 К с использованием методики, описанной в разделе 2.3. Аналогично, на рис. 22 показана восстановленная функция распределения частиц по радиусам для АТ= 20 К. Хорошо виден провал в области малых радиусов, обусловленный отсутствием в осевшей фракции частиц с малыми радиусами. Напомним, что ссдиментационная кривая представляет собой массовую долю осевшего вещества в зависимости от времени. Взвешенная фракция не дает вклада в седимептациоиную кривую. Таким образом, при использовании методики раздела 2.3, наличие узкого провала в функции распределения в области малых эквивалентных радиусов может служить указанием на то, что седиментационная картина искажается конвективными движениями среды. Работа опубликована в (Parchevsky, 2001а). Пусть нам дано N измерений функции J{x) на отрезке [а, Ь], содержащих экспериментальные ошибки. Пусть требуется найти первую (в общем случае k-ю) производную функции fix). Эту задачу можно решать различными методами. Однако, следует сразу отметить, что определение производной дх Лг- Л не может быть использовано для ее аппроксимации. Дело в том, что в эксперименте измеряется функция у(х) j{x) + с, а не/[д:). Здесь е- шум. Второе слагаемое вообще не имеет производной в классическом смысле. Если мы все же попытаемся аппроксимировать производную как Ayf&x, то получим пилообразную знакопеременную функцию большой амплитуды и реальная производная f{x) «утонет» в шуме, генерируемом вторым слагаемым. Восстановление производной из зашумленных данных является некорректной обратной задачей и, так или иначе, требует для своего решения специальных регуляризирующих. процедур. Существует несколько подходов к решению этой проблемы. Можно аппроксимировать J[x) какой-нибудь достаточно гладкой аналитической функцией, чтобы сгладить «шероховатости» от инструментальных ошибок, после чего производная может быть вычислена аналитически. Выбор вида аппроксимирующей функции очень важен. С одной стороны, аппроксимирующая функция должна содержать достаточно свободных параметров, чтобы иметь возможность аппроксимировать сложные функции. С другой стороны, она должна равномерно приближать искомую функцию па всем интервале. Эти требования накладывают достаточно жесткие условия на вид аппроксимирующей функции. Например, полиномиальная регрессия не подходит на роль аппроксимирующей функции, так как она не приближает искомую функцию равномерно. В качестве неплохой альтернативы может быть выбран кубический сплайн, который удовлетворяет обоим требованиям.Моделирование рассеяния акустических волн солнечным пятном
Определение аппаратной функции телескопа и восстановление изображений солнечных протуберанцев
Численное моделирование седиментации полидисперсной суспензии в присутствии двумерной сжимаемой конвекции
Восстановление скорости роста и функции плотности вероятности из зашумленных экспериментальных данных
Похожие диссертации на Применение численного моделирования и регуляризирующих алгоритмов для изучения конвекции и акустики в задачах физики Солнца и экологии
