Акустический журнал, 2022, T. 68, № 1, стр. 68-82

ВВЕДЕНИЕ

Дистанционная диагностика неоднородностей природных сред представляет значительный практический интерес, например, для оценки параметров грунта при строительстве и оценке устойчивости инженерных сооружений, обнаружения зарождающихся карстовых провалов, оползневых процессов и т.п. Значительный интерес представляет анализ характеристик приповерхностных слоев до глубин в несколько десятков метров [1], в частности, поиск таких потенциально опасных структур как карст. Существует множество методов дистанционного зондирования природных сред (электрическая и магнитная разведка, георадар, сейсмоакустическое зондирование). Примеры использования с обсуждением достоинств и недостатков (ограничений) различных геофизических методов исследования и диагностики неоднородностей представлены в обзорах [2, 3]. Отмечается, что сейсмоакустические методы характеризуются большей глубиной проникновения, а также меньшей чувствительностью к мелкомасштабным неоднородностям (загрязнению) вблизи поверхности и к наличию сильно искажающих электрические (электромагнитные) характеристики среды водных слоев. Примеры сейсмоакустических исследований на малых глубинах, где в качестве пробных волн используются объемные волны, для которых применяются стандартные для сейсмической разведки методы анализа отраженных и преломленных волн [4, 5], хорошо известны (например, [6–8]). Физические основы инженерной сейсмоакустики представлены, например, в [9–11].

Наряду с объемными интерес представляют поверхностные волны, характеристики которых также отражают особенности строения пористой среды. Наличие неоднородности среды по глубине приводит к появлению дисперсии поверхностной волны, что позволяет производить дистанционную диагностику слоистой структуры [12–14]. Адекватность модели в виде набора плоских параллельных слоев обусловлена накоплением осадочных отложений в поле тяжести Земли, что приводит к появлению подобных плоскослоистых структур. Теория распространения волн в слоистых структурах имеет долгую историю (например, [15, 16]), а обзор современного состояния проблемы представлен, например, в недавних монографиях [17, 18]. История и современное состояние методов диагностики на поверхностных волнах Рэлея и Стоунли представлена в [19, 20]. В последние примерно 10–20 лет активно развиваются методы, использующие природные источники шума. Впервые идея была высказана в пионерской работе Аки 1957 г. [21], но получила развитие в начале 21 века [22–24] с появлением совершенных цифровых систем регистрации, имеющих точную синхронизацию по времени. О растущем интересе к этой области исследований говорят и недавние публикации в Акустическом журнале [25, 26].

В представленной работе рассматриваются приложения инженерной сейсмоакустики с использованием портативного когерентного источника [27], что развивает идеи когерентной малоглубинной сейсмоакустической диагностики [8]. Достоинством методов диагностики, опирающихся на анализ поверхностных волн, является то, что в этом случае рассматривается распространение напрямую от источника к приемнику [28–30], что исключает зависимость амплитуды и фазы от угла падения волны на границу. Использование векторного приема позволяет определить параметр эллиптичности волны Рэлея (отношение амплитуд проекций горизонтального и вертикального смещения) и в результате получить оценку распределения коэффициента Пуассона по глубине [30]. Информация о величине коэффициента Пуассона важна для определения характера связей внутри гранулированной среды (см. модельные представления в [31]). Несмотря на ряд успешных примеров дистанционной диагностики с использованием векторного приема волны Рэлея [30, 32], остался до конца невыясненным вопрос о точности оценки параметров слоистой среды при использовании векторного приема и анализе эллиптичности волны Рэлея. Ответу на этот вопрос посвящена настоящая статья.

Точность оценки параметров определяется отношением амплитуды полезного или зондирующего сигнала к уровню шума [33]. Отметим, что в одной из первых работ по использованию многоканальных приемников для анализа дисперсии поверхностной волны Рэлея (например, [14]) вопрос о присутствии шума в данных не обсуждается. Тем самым предполагается, что амплитуда сигнала может быть задана произвольной, шум является несущественным фактором, и также не принимается во внимание возможная нелинейность среды в области возбуждения. Авторы настоящей статьи полагают важным учет отношения сигнал/шум в задаче оценки параметров слоистой среды методами SASW/MASW – аббревиатуры отвечают сокращениям: Spectral Analysis of Surface Waves/Multichannel Analysis of Surface Waves. При этом уровень шума целесообразно задавать, опираясь на результаты измерений. Отношение сигнал/шум определяется амплитудой зондирующего сигнала. В когерентной сейсмоакустике амплитуда задается малой настолько, чтобы обеспечить долговременную стабильность излучения [27]. В этом случае обеспечивается линейность режима возбуждения волн и отсутствие пластических деформаций в области возбуждения. Достижение необходимого соотношения сигнал/шум возможно за счет когерентного накопления многих посылок зондирующих сигналов. Поэтому еще одной целью настоящей работы является определение требований к характеристикам излучения для достижения требуемой точности реконструкции параметров среды.

Анализ задачи в самом общем случае большого смысла не имеет из-за многообразия реализуемых в природе конфигураций слоистых сред. На взгляд авторов, более правильным будет использование моделей, отвечающих каким-то типичным ситуациям. В практически значимых случаях анализа глубин до первых десятков метров в качестве модели среды или ее первого приближения достаточно рассмотреть модель из одного-двух слоев, лежащих поверх полупространства. Такая модель в условиях, типичных для Восточно-Европейской платформы, отвечает верхнему слою почвы, лежащему на суглинках или супеси поверх глинистых отложений. Например, для условий Нижегородской области консолидированные осадочные породы типа карбонатных или песчаников размещаются на глубине порядка $50$ метров [1], что позволяет упростить модель среды на меньших глубинах. Глубинный разрез в условиях сейсмоакустического полигона ИПФ РАН “Безводное” был получен методом фазового межскважинного профилирования на SH-волнах [28] и подтвержден результатами бурения с отбором керна. В этом случае было наглядно показано присутствие слоя почвы толщиной до 0.5 м и затем слоя супеси и суглинка с переходом к глине. При этом контраст по скоростным характеристикам слоев составил около 10–30% процентов. Таким образом, точность оценки параметров должна быть не хуже указанной величины для разрешения особенностей структуры поверхностных отложений и эффективного использования дистанционных методов диагностики сред в инженерной сейсмике.

Перечислим основные цели работы:

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

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

Работа организована следующим образом. В первой части дано краткое описание схемы измерений, используемой в методах SASW/MASW, а также способ представления данных. Во второй части приведены данные по измеренным спектральным амплитудам природного шума, а также рассмотрен вопрос о точности измерения параметров слоистой среды. В третьей части на примере модели в виде слоя, лежащего на полупространстве, доказывается возможность измерения параметров слоистой структуры с требуемой для целей диагностики точностью. Здесь же показывается, что привлечение параметра эллиптичности обеспечивает оценку коэффициента Пуассона слоев с удовлетворительной для целей диагностики точностью. В третьей части работы приведен пример структуры из трех слоев, где промежуточный слой имеет отличия только в коэффициенте Пуассона по сравнению с вмещающими его слоями. Такая модель отвечает практически важному случаю наличия насыщенного жидкостью пористого слоя между непроницаемыми средами сверху и снизу. Здесь показано, что эффективность диагностики насыщения повышается за счет использования частотной зависимости параметра эллиптичности. Все приведенные примеры демонстрируют практическую реализуемость развиваемых авторами методов дистанционной диагностики природных сред.

СХЕМА ИЗМЕРЕНИЙ И СПОСОБ ПРЕДСТАВЛЕНИЯ ДАННЫХ

На рис. 1 представлена типичная схема сейсмоакустических измерений, целью которых является определение параметров слоев среды. Вибрационный источник подробно описан в работе [27], его фотография приведена справа на рис. 1. На этой фотографии видны два акселерометра, с помощью которых осуществляется определение силы, действующей на грунт. Эта процедура подробно описана в [27] и в дальнейшем не упоминается: важно, что имеется однозначная и известная связь между амплитудой электрического сигнала и силой, возбуждающей сейсмические волны. В качестве зондирующего используется электрический сигнал с линейной частотной модуляцией (ЛЧМ) [34], который для подавления переходных процессов при включении и выключении вибратора подвергается процедуре взвешивания. В начальный и конечный моменты времени излучения сам электрический сигнал и его производная по времени равны нулю. В качестве весовой функции в наших экспериментах [27–30, 35] используется косинусное окно [36].

Отличие между многоканальным (MASW) [13] и двухканальным (SASW) [12] методами спектрального анализа поверхностных волн фактически сводятся к повышению отношения сигнал/шум за счет когерентной обработки данных с нескольких приемников (выигрыш от использования линейной антенны из геофонов). Поэтому при дальнейшем рассмотрении мы опираемся на простейшую реализацию с использованием двух геофонов [12]. В этом случае оценка скорости волны Рэлея ${{V}_{R}}$ выполняется по разности фаз $\Delta \varphi $ сигналов двух геофонов, установленных на линии, соединяющей их с источником. Параметр эллиптичности равен отношению спектральных амплитуд соответствующих ненулевых проекций вектора колебательной скорости (смещения, ускорения) дневной поверхности:

где $d$ – расстояние между геофонами, $\omega $ – круговая частота, ${{G}_{{12}}}(\omega )$ – взаимный спектр [36] для двух приемников, регистрирующих одну из проекций вектора колебательной скорости (смещения, ускорения), $G_{{11}}^{{(rz)}}(\omega )$ и $G_{{22}}^{{(rz)}}(\omega )$ – взаимные спектры каждого из приемников для двух проекций вектора колебательной скорости (смещения, ускорения), $G_{{11}}^{{(zz)}}(\omega )$ и $G_{{22}}^{{(zz)}}(\omega )$ – спектры мощности проекции $z$ двух приемников. Величина $\Delta \varphi (\omega )$ определяется как фаза взаимного спектра двух приемников ${{G}_{{12}}}(\omega )$ [12]. При этом предполагается, что приемник 1 расположен ближе к источнику, чем приемник 2. Для оценки скорости волны Рэлея (1) неважно, какие проекции вектора колебательной скорости использованы, важным является определение разности фаз между одинаковыми проекциями, что отражают верхние индексы в выражении для взаимного спектра. Оценка параметра эллиптичности (2) выполняется по отношению спектральных амплитуд соответствующих проекций двух векторных приемников. При распространении волны Рэлея в плоскослоистой среде от нуля отличаются две проекции: ${{v}_{r}}$ – горизонтальная вдоль линии, соединяющей приемники с точкой приложения силы, и ${{v}_{z}}$ – вертикальная.

Точность измерений определяется двумя факторами. Во-первых, очевидным влиянием посторонних шумов. Во-вторых, особенностями излучения пробных сейсмоакустических волн. Оценку влияния второго фактора лучше всего иллюстрирует результат численного моделирования для однородного полупространства. В этом случае волна Рэлея не имеет дисперсии. Точность измерения скорости (1) грубо определяет минимально обнаруживаемый контраст параметров среды. Анализ максимально простой модели среды позволяет определить методические проблемы, возникающие при использовании спектрального метода анализа поверхностной волны и когерентных источников сейсмоакустического излучения. На рис. 2 внизу представлен результат численного моделирования оценки скорости волны Рэлея (1) при использовании линейно-модулированного по частоте сигнала силы в полосе частот 10–100 Гц и расстоянии между двумя геофонами $d = 1$ м. Сверху представлен результат оценки (2) отношения проекций вектора колебательной скорости на дневной поверхности (эллиптичности волны Рэлея).

Искажения данных на низких частотах и частотах вблизи $100$ Гц обусловлены наличием весового косинусного окна. Точное значение скорости волны Рэлея на границе однородного упругого полупространства равно ${{V}_{R}} = 0.9325{{V}_{S}}$, где скорость волны сдвига была задана равной ${{V}_{S}} = 150$ м/с. Конкретное численное значение отношения ${{{{V}_{R}}} \mathord{\left/ {\vphantom {{{{V}_{R}}} {{{V}_{S}}}}} \right. \kern-0em} {{{V}_{S}}}}$ имеет слабую зависимость от коэффициента Пуассона [37] и отвечает ${{{{V}_{S}}} \mathord{\left/ {\vphantom {{{{V}_{S}}} {{{V}_{P}}}}} \right. \kern-0em} {{{V}_{P}}}} = {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$, где ${{V}_{P}}$ – скорость продольной волны в упругом полупространстве. Заданные в численной модели значения скоростей объемных волн ${{V}_{S}}$ и ${{V}_{P}}$ отвечают многолетним измерениям в условиях полигона “Безводное” ИПФ РАН [27, 35].

Ожидаемое значение скорости волны Рэлея равно $139.875$ м/с и показано на нижнем графике рис. 2 штриховой линией. На верхнем графике рис. 2 штриховой линией показано теоретическое значение отношения амплитуд проекций $E \approx 0.639$, отвечающее ${{{{V}_{S}}} \mathord{\left/ {\vphantom {{{{V}_{S}}} {{{V}_{P}}}}} \right. \kern-0em} {{{V}_{P}}}} = {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$ или величине коэффициента Пуассона $\nu = {1 \mathord{\left/ {\vphantom {1 3}} \right. \kern-0em} 3}$. Заштрихованные области на графиках показывают отклонения от теоретических значений на ±0.1%. Скорость волны Рэлея (1) в полосе частот 20–90 Гц, где весовое окно вносит минимальные возмущения в спектр, лежит в пределах заштрихованной области. В силу линейной связи между скоростью волны Рэлея и скоростью сдвиговой волны последняя определяется с той же относительной точностью около 0.1%. Параметр эллиптичности (2) не связан с измерением фазы, и поэтому определяется с более высокой точностью. Таким образом, в отсутствие шума особенности возбуждения волн проявляются в области частот, где заметно отклонение весовой функции от единицы.

ХАРАКТЕРНЫЙ УРОВЕНЬ ШУМА, ПОГРЕШНОСТЬ ОЦЕНКИ ПАРАМЕТРОВ

Для оценки влияния шума на точность измерений необходимо располагать информацией о характерных значениях спектральной плотности амплитуды природного шума. Эти амплитуды предполагается в дальнейшем использовать в численном моделировании. На рис. 3 представлены измеренные в различное время и с использованием различных средств спектральные плотности амплитуды шума вертикальной проекции вектора ускорения дневной поверхности в районе полигона “Безводное” ИПФ РАН. Тонкие серая и черная линии отвечают измерениям 2006 года на полигоне, где проводились эксперименты [27, 35]. Место проведения экспериментов удалено от автомобильных дорог и иных источников индустриального шума и представляет собой открытый участок в лесу. Для измерений использовались акселерометры и усилители заряда на операционных усилителях с низким уровнем теплового шума. Кривая для теплового шума также приведена на графике серой линией с длинным штрихом, чтобы иметь возможность судить о достоверности измерений природного сейсмического шума. Хорошо видно, что этот шум в области частот выше нескольких герц обусловлен ветром, который возбуждает сейсмические колебания, например, через корневую систему колеблющихся деревьев.

Жирные сплошная и штриховая линии отвечают измерениям с использованием векторных сейсмических станций Zet048C, содержащих три геофона с рабочим диапазоном частот 0.05–50 Гц (более подробная информация представлена на сайте производителя оборудования www.zetlab.ru). Уровень собственного шума станций близок к минимальному наблюдаемому уровню, который показан на рис. 3 линией с символами и в литературе имеет аббревиатуру NLNM [38]. Данные рис. 3 отвечают шумовой компоненте вертикального ускорения дневной поверхности в двух записях. Сплошная линия, обозначенная в подписи на рис. 3 как Zetlab 1, отвечает ночной записи ноября 2018 г. при установке станции на бетонное основание в подвале капитального строения, удаленного от антропогенных источников. Жирная темно-серая линия, обозначенная в подписи на рис. 3 как Zetlab 2, отвечает записи длительностью около месяца, которая была выполнена на полигоне “Безводное”.

Измерениям отвечает низкий уровень низкочастотного шума с хорошо различимым максимумом в окрестности $0.2$ Гц, который связан с генерацией штормовых микросейсм [39]. В момент проведения измерений наблюдались порывы северо-восточного ветра, что проявилось в увеличении уровня шума для соответствующих проекций. Отметим, что данные, полученные с помощью сейсмостанции Zet048C, согласуются с данными, полученными с помощью акселерометров (первые три линии), и позволяют судить о возможных вариациях природного шума в полосе частот от единиц до сотен герц.

В целом шум в области частот выше $1$ Гц, предположительно, в наибольшей степени связан с ветром, поскольку расстояние до ближайшей автомобильной трассы составляло примерно 20 км, а иных близлежащих антропогенных источников не имелось при проведении всех измерений, данные которых представлены на рис. 3. Исследование особенностей формирования сейсмического шума в условиях Восточно-Европейской платформы представлено в [40], где авторы указывают на антропогенный характер сейсмического шума в области частот выше $1$ Гц, а также значительную корреляцию шума в этой области частот с силой ветра. Линия с символами отвечает минимальному среднему уровню шума по мировой сети сейсмических станций [38]. Эта кривая может рассматриваться в качестве минимально возможной (наблюдаемой) спектральной амплитуды природного шума. Подобный уровень шума возможно наблюдать в малонаселенных местах (например, [41]), к каковым, очевидно, не относится Восточно-Европейская платформа. Выполненные нами измерения указывают на значительное превышение уровня шума по сравнению с уровнем NLNM. Отметим, что и амплитуды шума, измеренные авторами [40] в Воронежской области, также существенно превышают уровень NLNM. Горизонтальные штриховые линии со стрелками показывают характерные минимальные и максимальные уровни шума, а заштрихованная область отвечает диапазону изменений.

Диапазон частот в сейсмоакустических исследованиях определяется, во-первых, характеристиками источника, во-вторых, полосой рабочих частот регистратора и приемника. Значительный интерес представляют вибрационные источники электромагнитного типа [6, 27]. В этом случае источник характеризуется высокой стабильностью параметров излучаемых сейсмических волн [27], что обеспечивает долговременное накопление и повышение отношения сигнал/помеха. Подобные портативные вибрационные источники имеют нижнюю рабочую частоту порядка 10 Гц. Многоканальные сейсмические станции имеют рабочий диапазон частот от нескольких герц до первых килогерц. В качестве примера можно привести сейсмостанцию Лакколит Х-М2 с полосой рабочих частот 5–8000 Гц. Распространенные в инженерной сейсмике геофоны, например, GS–20DX имеют полосу рабочих частот 10–250 Гц при незначительных вариациях верхней и нижней рабочих частот для других типов геофонов. На рис. 3 штрих-пунктирной линией показан уровень электронного шума сейсмостанции Лакколит Х-М2. Нетрудно видеть, что на частотах ниже 100 Гц этот шум заметно ниже природного сейсмоакустического шума, и электронный шум сейсмостанции не приводит к ограничению точности измерения параметров слоистых сред методами SASW/MASW.

Расширение диапазона рабочих частот в область высоких частот большого смысла не имеет, поскольку на 200 Гц длина поверхностной волны Рэлея для условий полигона “Безводное” составляет порядка 0.5 метра. Самый верхний слой оказывается сильно неоднородным, что приводит к значительному поглощению волн из-за рассеяния (данные приведены в работе [35]). Использование более низких частот неизбежно приводит к увеличению массы и габаритов геофонов и вибрационных источников, что неудобно для оперативного обследования площадок под строительство, локализации нарушений в структуре среды. Таким образом, несложные оценки показывают, что при характерной величине скорости сдвиговой волны ~200 м/с [27, 35], отвечающей типичным значениям для приповерхностных слоев в средней полосе России (например, данные в [30]), частотный диапазон 5–200 Гц регистрирующей аппаратуры в полной мере обеспечивает исследования природных сред на глубинах от десятка сантиметров до примерно 20 метров с использованием поверхностной волны Рэлея в качестве зондирующей. Таким образом, диапазон рабочих частот составляет 10–200 Гц с возможными небольшими вариациями верхней и нижней границ. Этот диапазон частот и вариации уровня природного сейсмического шума показаны на рис. 3 заштрихованной светло-серым цветом областью.

На рис. 4 представлена оценка функции распределения шумового сигнала вертикальной проекции вектора колебательной скорости в полосе частот выше 1 Гц. Число интервалов группировки ${{k}_{g}}$ при оценке плотности вероятности выбиралось в соответствии с правилом [36]: ${{k}_{g}} = 1.87{{(N - 1)}^{{2/5}}} \approx 717$, где $N = 2\,880\,001$ число временных отсчетов анализируемой записи. Правая ось ординат на рис. 4 показывает, что число “попаданий” в каждую ячейку гистограммы много больше десяти (общепринятое требование при построении гистограмм не менее 10 [42]). Серой линией показана аппроксимация функцией распределения нормального процесса с нулевым средним и стандартным отклонением $0.18$ мкм/с. Нетрудно видеть, что шум имеет статистику, очень близкую к нормальному распределению. Критерий проверки ${{\chi }^{2}}$ нормальности распределения (см., например, [36, 42]) указывает на не более, чем 5% отклонение от нормального закона распределения в пределах $\left| {{{v}_{z}}} \right|\,\,\leqslant \,\,0.8$ мкм/с. Отклонение от нормального распределения при больших амплитудах скорости, предположительно, связано с нестационарным характером шума (например, из-за наблюдавшихся в процессе измерений порывов ветра). Возможность использования нормального распределения для описания сейсмоакустического шума используется при оценке погрешности определения параметров слоистой среды методом спектрального анализа поверхностных волн.

Результат расчета на рис. 2 отвечает отсутствию шума, т.е. бесконечной величине отношения сигнал/шум, а отклонения величин ${{V}_{R}}$ и $E$ от теоретических значений обусловлены особенностями возбуждения пробных сейсмоакустических волн. Наличие естественных шумов (рис. 3), очевидно, должно приводить к деградации точности оценки $\delta {{V}_{R}}$. Линейно-модулированный по частоте сигнал (ЛЧМ) имеет постоянную спектральную плотность мощности в полосе частот излучения при условии большой величины базы сигнала $\mathfrak{m}$, равной произведению девиации частоты на длительность излучения [34]. Именно такие варианты использования ЛЧМ зондирующих сигналов применяются в когерентной сейсмоакустике, поскольку обеспечивают эффективное использование процедур свертки принятого сигнала с опорным [27, 34]. Численному примеру на рис. 2 отвечает $\mathfrak{m} = $ 738, т.е. $\mathfrak{m} \gg 1$. Используя выражения из статьи [43], можно получить следующую инженерную оценку спектральной плотности амплитуды вертикальной проекции вектора колебательной скорости на дневной поверхности:

где числовой коэффициент 0.22 отвечает типичному отношению скоростей объемных волн ${{{{V}_{S}}} \mathord{\left/ {\vphantom {{{{V}_{S}}} {{{V}_{P}}}}} \right. \kern-0em} {{{V}_{P}}}} = {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$, $\hat {F}(\omega )$ – спектральная плотность амплитуды вертикальной силы, возбуждающей сейсмические колебания, $\rho $ – плотность упругого полупространства, ${{H}_{0}}({{k}_{R}}r)$ – функция Ханкеля нулевого порядка, ${{k}_{R}} = {\omega \mathord{\left/ {\vphantom {\omega {{{V}_{R}}}}} \right. \kern-0em} {{{V}_{R}}}} = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } {{{\lambda }_{R}}}}} \right. \kern-0em} {{{\lambda }_{R}}}}$, ${{\lambda }_{R}}$ – длина волны Рэлея, $r$ – расстояние от точки приложения силы до приемника. В волновой зоне источника при условии $r\,\, \geqslant \,\,{{\lambda }_{R}}$ функцию Ханкеля с высокой точностью можно аппроксимировать асимптотическим выражением: $\left| {{{H}_{0}}({{k}_{R}}r)} \right| \approx \frac{1}{\pi }\sqrt {\frac{{{{\lambda }_{R}}}}{r}} $. В результате амплитуда колебательной скорости оказывается пропорциональной частоте в степени 3/2.

Природный сейсмоакустический шум характеризуется примерным постоянством спектральной плотности ускорения в полосе частот от единиц до сотен герц (рис. 3). Поэтому отношение сигнал/шум пропорционально частоте в степени 5/2. Это указывает на наличие сложностей в измерении дисперсии волны Рэлея в области низких частот. Еще раз отметим, что увеличение амплитуды силы вибрационного источника имеет естественные ограничения, обусловленные появлением пластических деформаций в области возбуждения и нестабильностью источника. При характерной для стабильной работы вибрационного источника величине силы $10$ Н [27] и полосе излучения 10–100 Гц при условии $\mathfrak{m} \gg 1$ спектральная плотность амплитуды силы равна $\hat {F}(\omega ) \approx 1$ Н/$\sqrt {{\text{Гц}}} $. Несложные оценки по формуле (3) показывают, что при такой спектральной амплитуде источника амплитуда ${{\hat {v}}_{z}}(\omega )$ в области частот порядка $10$ Гц оказывается сопоставимой со спектральной амплитудой природного сейсмоакустического шума.

На рис. 5 показан результат расчета частотной зависимости спектральной плотности ускорения ${{\hat {a}}_{z}}(\omega ) = \omega {{\hat {v}}_{z}}(\omega )$ на расстоянии 15 метров от источника. Приведение данных для ускорения обусловлено примерным постоянством спектральной плотности амплитуды шумового ускорения дневной поверхности (рис. 3). Черная линия отвечает излучению в полосе 10–100 Гц. Серые линии показывают очевидный выигрыш в отношении сигнал/шум при сокращении полосы излучения до 1 Гц ($\mathfrak{m} \approx 8$), т.е. при десятикратном увеличении спектральной плотности амплитуды источника. Отношение амплитуд узкополосных сигналов приблизительно отвечает отмеченной выше зависимости $ \propto {{\omega }^{{5/2}}}$. Уровень шума на рис. 5 отвечает зарегистрированному минимуму при нулевой скорости ветра (рис. 3).

В случае аддитивного шума относительные вариации амплитуды и фазы сигнала одинаковы. Следовательно, возмущения фазы составят порядка $\delta \varphi = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } {{{N}_{s}}}}} \right. \kern-0em} {{{N}_{s}}}}$, где ${{N}_{s}}$ – отношение амплитуды сигнала к амплитуде шума. Например, при использовании широкополосного излучения (жирная черная линия на рис. 5) на частоте 20 Гц отношение амплитуды возбуждаемой волны Рэлея к амплитуде шума составляет порядка ${{N}_{s}} = 10$ и $\delta \varphi \approx 0.2\pi $. Из выражения (1) несложно получить связь $\delta \varphi $ с погрешностью оценки скорости волны Рэлея $\delta {{V}_{R}}$:

Отсюда следует, что расстояние $d$ целесообразно задавать максимально возможным при сохранении однозначности фазы: $d \sim {{\lambda }_{R}}$. В этом случае измерения для условий рис. 5 дадут погрешность оценки порядка отношения сигнал/шум: ${{\delta {{V}_{R}}} \mathord{\left/ {\vphantom {{\delta {{V}_{R}}} {{{V}_{R}}}}} \right. \kern-0em} {{{V}_{R}}}} \sim {1 \mathord{\left/ {\vphantom {1 {{{N}_{s}}}}} \right. \kern-0em} {{{N}_{s}}}}$. Погрешность измерения отношения $\delta E = {1 \mathord{\left/ {\vphantom {1 {{{N}_{s}}}}} \right. \kern-0em} {{{N}_{s}}}}$, что следует из (1). Погрешность оценки коэффициента Пуассона зависит от производной ${{dE} \mathord{\left/ {\vphantom {{dE} {d\nu }}} \right. \kern-0em} {d\nu }}$. Зависимости $E$ и ${{dE} \mathord{\left/ {\vphantom {{dE} {d\nu }}} \right. \kern-0em} {d\nu }}$ от коэффициента Пуассона приведены на рис. 6. В представляющей интерес области изменения $0\,\,\leqslant \,\,\nu < {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$ модуль производной $\left| {{{dE} \mathord{\left/ {\vphantom {{dE} {d\nu }}} \right. \kern-0em} {d\nu }}} \right| \approx {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$. Следовательно, погрешность оценки коэффициента Пуассона имеет порядок:

т.е. погрешность оценки коэффициента Пуассона имеет тот же порядок, что и относительная погрешность оценки скорости волны Рэлея: $\delta \nu \sim {{\delta {{V}_{R}}} \mathord{\left/ {\vphantom {{\delta {{V}_{R}}} {{{V}_{R}}}}} \right. \kern-0em} {{{V}_{R}}}}$.

Увеличение точности возможно за счет роста ${{N}_{s}}$ при использовании когерентного накопления импульсных откликов для каждого сеанса излучения. Отношение сигнал/шум растет пропорционально квадратному корню из числа ${{N}_{p}}$ сеансов излучения ЛЧМ сигнала: ${{N}_{s}} \propto \sqrt {{{N}_{p}}} $. Для достижения точности порядка 1% (в десять раз шире заштрихованной области на рис. 2) необходимо использовать порядка 100 сеансов излучения или примерно 14 мин измерений. Отметим, что использование многих приемников (многоканальная реализация спектрального анализа поверхностных волн – MASW) аналогично когерентному усреднению временных реализаций. Например, при использовании вместо двух приемников ${{N}_{g}}$ геофонов выигрыш в отношении сигнал/шум составит $\sqrt {{{{{N}_{g}}} \mathord{\left/ {\vphantom {{{{N}_{g}}} 2}} \right. \kern-0em} 2}} $ раз. Данные рис. 5 явно указывают на еще одну возможность увеличения отношения сигнал/шум – сужение полосы частот ЛЧМ сигнала. Это следует из обратной пропорциональности спектральной плотности мощности полосе излучаемых частот при большой базе сигнала. При этом в каждом конкретном случае необходимо убедиться в возможности подавлении вкладов волн, отличных от волны Рэлея [35], для сокращенной полосы излучения.

ТОЧНОСТЬ ОЦЕНКИ ПАРАМЕТРОВ ЗАГЛУБЛЕННОГО СЛОЯ

Рассмотрим простую модельную среду, состоящую из упругого слоя мощности $D = 1$ м с параметрами ${{V}_{{1S}}} = 150$ м/с и ${{V}_{{1P}}} = 300$ м/с, лежащего поверх полупространства с ${{V}_{{2S}}} = 160$ м/с и ${{V}_{{2P}}} = 320$ м/с. Такая структура приблизительно отвечает условиям полигона “Безводное” ИПФ РАН до глубин порядка $5$ метров [28, 27 ] и является достаточно типичной для Нижегородской области [1, 30]. Контраст изменения скорости сдвиговой волны сознательно задан малым (примерно 7%), чтобы показать, с одной стороны, сложность диагностики такой неоднородности, а с другой, – возможность определения параметров среды средствами когерентной сейсмоакустики. Если рассматривать простейшую реализацию метода SASW/MASW, которая не включает анализ проекций вектора смещения дневной поверхности, то величины плотности $\rho $ и продольной скорости ${{V}_{P}}$ не подлежат определению, как не влияющие на дисперсию волны Рэлея [14, 20]. Светло-серым цветом на рис. 7 заштрихована область вариаций измеряемых величин при изменении скорости продольной волны в пределах $270\,\,\leqslant \,\,{{V}_{{2P}}}\,\,\leqslant \,\,370$ м/с. Нетрудно видеть, что вариации дисперсионной зависимости существенно меньше, чем изменения параметра эллиптичности. Таким образом, при анализе только дисперсионной зависимости задача может состоять только в определении скоростей волны сдвига в слоях и глубины раздела слоев при заданных отношениях ${{{{V}_{{1S}}}} \mathord{\left/ {\vphantom {{{{V}_{{1S}}}} {{{V}_{{1P}}}}}} \right. \kern-0em} {{{V}_{{1P}}}}}$ и ${{{{V}_{{2S}}}} \mathord{\left/ {\vphantom {{{{V}_{{2S}}}} {{{V}_{{2P}}}}}} \right. \kern-0em} {{{V}_{{2P}}}}}$ (см. также обсуждение в [20]).

На рис. 7 сплошной черной линией показан результат расчета по формуле (1) в отсутствие шума. Нетрудно видеть, что погрешность определения скорости волны Рэлея в отсутствие шума примерна такая же, как и для однородного полупространства. Заметные отклонения, как и для рис. 2, связаны с наличием весового окна, исключающего переходные процессы при работе излучателя [27]. Таким образом, все результаты, которые получены для простейшей модели однородного полупространства, имеют место и в более сложной модели двух слоев малого контраста скоростных параметров. Серой линией показан результат расчета по формуле (1) при наличии шума. Отношение сигнал/шум отвечает данным рис. 5 (моделирование измерений в тихом месте при слабом ветре).

На рис. 7 вверху показаны частотные зависимости отношения амплитуд проекций вектора колебательной скорости (2). Рукав в виде областей светло-серого цвета показывает вариации дисперсионной зависимости и $E$ при изменении скорости продольной волны ${{V}_{{2P}}}$ на $ \pm 50$ м/с (т.е. при вариации на 16%). Относительные изменения дисперсионной зависимости не превышают 1.4%, т.е. существенно меньше заданного контраста скорости сдвиговой волны. В то же время относительные изменения величины $E$ составляют величину около 7%, что имеет приблизительно тот же порядок. Предложение использовать анализ частотной зависимости $E$ было высказано и проверено экспериментально в работах [30, 32]. Ниже рассмотрен опущенный ранее вопрос о точности оценки коэффициента Пуассона при анализе частотной зависимости величины эллиптичности (2).

Задача определения параметров методами спектрального анализа поверхностных волн рассмотрена в [14] и сводится к модификации градиентного метода Ньютона. Поскольку остаточная погрешность определяется уровнем шума, то при условии нормального распределения шума (рис. 4) точность определения параметров равна [33]:

где $\Delta {{p}_{j}}$ представляет собой стандартное отклонение величины ${{p}_{j}}$, ${{\epsilon }_{1}}$ – взвешенная остаточная дисперсия между измеренными величинами скорости волны Рэлея и их модельным представлением, ${\mathbf{p}}$ – вектор параметров модели. При анализе только дисперсионной зависимости величина ${{\epsilon }_{1}}$ имеет размерность квадрата скорости, вектор параметров ${\mathbf{p}} = \{ {{V}_{{1S}}},{{V}_{{2S}}},D\} $ представляет собой скорости сдвиговых волн и глубину залегания границы слоев (в общем случае, набор скоростей и глубин залегания границ слоистой структуры).

Величина $H_{{ij}}^{{(1)}}$ представляет собой матрицу Гесса:

где

Здесь ${{y}_{m}}$ – вычисленная скорость волны Рэлея ${{V}_{R}}$ на частоте ${{f}_{m}}$ из набора частот ${{f}_{1}} - {{f}_{M}}$, на которых выполнены измерения этой скорости. Величина ${{w}_{m}}$ представляет собой безразмерный весовой множитель, который используется для учета возможной неравномерной погрешности измерений для каждого отсчета, что может быть связано с зависимостью ${{N}_{s}}$ от частоты (рис. 5). В качестве веса данных часто используется величина, обратная дисперсии. При этом сумма весовых коэффициентов должна быть равна единице. В этом случае величины ${{w}_{m}}$ можно определить следующим образом:

Весовые множители (8) учитывают отношение мощности шума к мощности зондирующего сейсмоакустического сигнала. При вычислении величины остаточной дисперсии также используются коэффициенты ${{w}_{m}}$: ${{\epsilon }_{1}} = \sum\nolimits_{m = 1}^M {{{w}_{m}}} {{\left( {{{y}_{m}} - {{{\tilde {y}}}_{m}}} \right)}^{2}}$, где ${{\tilde {y}}_{m}}$ – измеренные величины скорости волны Рэлея.

Учет частотной зависимости отношения амплитуд проекций вектора колебательной скорости (2) также приводит к выражениям, аналогичным (6) и (7):

где $\Delta {{\nu }_{j}}$ представляет собой стандартное отклонение величины коэффициента Пуассона ${{\nu }_{j}}$ в слое с номером $j$, ${{\epsilon }_{2}}$ – взвешенная остаточная дисперсия между измеренными величинами $E$ и их модельными значениями (весовые коэффициенты определены выражением (8)). Поскольку отношение амплитуд проекций является безразмерной величиной, ${{\epsilon }_{2}}$ также безразмерная величина. Матричные элементы $H_{{ij}}^{{(2)}}$ равны: где

Здесь ${{y}_{m}}$ – вычисленное отношение $E$ на частоте ${{f}_{m}}$ из набора частот ${{f}_{1}} - {{f}_{M}}$, на которых выполнены измерения этой величины. Величина ${{w}_{m}}$ определена (8).

Величины ${{\partial {{y}_{m}}} \mathord{\left/ {\vphantom {{\partial {{y}_{m}}} {\partial {{p}_{j}}}}} \right. \kern-0em} {\partial {{p}_{j}}}}$, входящие в (7), и ${{\partial E} \mathord{\left/ {\vphantom {{\partial E} {\partial {{\nu }_{j}}}}} \right. \kern-0em} {\partial {{\nu }_{j}}}}$, входящие в (10), показаны на рис. 8. Максимум величины $\left| {{{\partial {{V}_{R}}} \mathord{\left/ {\vphantom {{\partial {{V}_{R}}} {\partial D}}} \right. \kern-0em} {\partial D}}} \right|$ достигается при приближенном равенстве ${{{{\lambda }_{R}}} \mathord{\left/ {\vphantom {{{{\lambda }_{R}}} 2}} \right. \kern-0em} 2} \approx D$, что согласуется с результатом работы [14]. Штриховыми линиями на рис. 8 показаны производные по скоростям продольных волн в слоях. Эти величины примерно на 1–2 порядка меньше производных по скоростям сдвиговых волн. Приведенные зависимости иллюстрируют известный факт (например, [20]) практической невозможности оценки скоростей продольных волн в стандартных приложениях SASW/MASW, использующих анализ только вертикальной проекции вектора колебательной скорости дневной поверхности.

Величины производных от параметра эллиптичности по коэффициенту Пуассона имеют порядок единицы, что ясно указывает на возможность использования зависимости $E(\omega )$ для оценки коэффициента Пуассона слоистой среды. При этом оценка глубины $D$ залегания границы на основе анализа частотной зависимости $E$ нецелесообразна из-за малой величины соответствующей производной ${{\partial E} \mathord{\left/ {\vphantom {{\partial E} {\partial D}}} \right. \kern-0em} {\partial D}}$ по сравнению с производной ${{\partial {{V}_{R}}} \mathord{\left/ {\vphantom {{\partial {{V}_{R}}} {\partial D}}} \right. \kern-0em} {\partial D}}$.

Погрешности определения параметров, вычисленные по формулам (6) и (9), составляют:

Погрешности определения скоростей ${{V}_{{1S}}}$ и ${{V}_{{2S}}}$ существенно меньше заданного контраста скоростей. Наихудшая погрешность наблюдается для оценки коэффициента Пуассона второго (заглубленного) слоя. Этой величине погрешности отвечает изменение скорости продольной волны $\Delta {{V}_{{2P}}} = \pm 50$ м/с, отмеченное серой штриховкой на верхнем графике (рис. 7).

Из представленных на рис. 8 зависимостей следует стратегия поиска параметров слоистой среды спектральным методом анализа поверхностных волн SASW/MASW при учете параметра эллиптичности. Положение границ и оценка скорости волн сдвига в слоях выполняется также, как и в ранее (например, [14]) при заданном из априорных соображений коэффициенте Пуассона [20]. Анализ частотной зависимости $E(\omega )$ позволяет получить уточненное значение коэффициента Пуассона для слоев на каждом шаге итерационной процедуры [14]. Из-за наличия слабой зависимости скорости волны Рэлея от коэффициента Пуассона полученная оценка приведет к появлению небольших поправок к оценкам скоростей волн сдвига. Последовательное использование итераций, очевидно, как и в случае стандартной реализации SASW/MASW, приведет к получению оценки обеих скоростей в слоях.

Возможность дистанционной оценки величины коэффициента Пуассона важна для определения структуры связей между зернами в гранулированных средах (модели сред описаны, например, в [31]). Величина $\nu \approx 0$ отвечает контакту Герца–Миндлина, при котором наряду с прижимающими зерна друг к другу силами действуют силы сухого трения, препятствующие их относительному сдвигу. Величина $\nu = {1 \mathord{\left/ {\vphantom {1 4}} \right. \kern-0em} 4}$ отвечает модели центральных сил между зернами, когда разрешено относительное проскальзывание зерен. Величина $\nu = {1 \mathord{\left/ {\vphantom {1 3}} \right. \kern-0em} 3}$ (отношение скоростей объемных волн ${{{{V}_{S}}} \mathord{\left/ {\vphantom {{{{V}_{S}}} {{{V}_{P}}}}} \right. \kern-0em} {{{V}_{P}}}} = {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$) отвечает типичному значению коэффициента Пуассона для большинства упругих материалов. Наконец, величина $\nu \approx {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}$ отвечает гранулированным средам, пористое пространство в которых заполнено жидкостью. Таким образом, приведенная выше погрешность оценки коэффициента Пуассона $\delta \nu = 0.12$ позволяет разрешить перечисленные основные варианты структуры связей в гранулированных средах.

ПРИМЕР МОДЕЛЬНОЙ СРЕДЫ, В КОТОРОМ ВАЖНА ОЦЕНКА КОЭФФИЦИЕНТА ПУАССОНА ДЛЯ ДИАГНОСТИКИ НАСЫЩЕНИЯ ПОР ЖИДКОСТЬЮ

В качестве еще одного примера рассмотрим модель среды, где верхний слой и упругое полупространство ограничивают полностью насыщенный водой пористый слой. Параметры модельной среды из трех упругих слоев зададим таким образом, чтобы скорости волны Рэлея для каждой из трех сред были одинаковыми. Покажем, что использование частотной зависимости $E(\omega )$ позволяет разрешить насыщенный водой слой, а при анализе только дисперсионной зависимости скорости волны Рэлея, что характерно для известных методов SASW/MASW, получается неверное решение.

Параметры модельных сред приведены в табл. 1. Модель 1 отвечает параметрам пористой среды с заполнением пор водой. Модель 2 отвечает среде с постоянным значением коэффициента Пуассона для всех слоев. При этом параметры заданы таким образом, что дисперсионные свойства волны Рэлея в моделях 1 и 2 практически неразличимы в ограниченном частотном диапазоне (в данном примере выше $10$ Гц). Поскольку всегда из-за особенностей используемых источников и регистрирующей аппаратуры (Приложение) присутствуют ограничения снизу на диапазон частот, в пределах которого производятся измерения, рассматриваемый нами пример не является каким-то частным случаем.

Величина $D$ отвечает глубине нижней границы слоя, прочерк для последнего означает полупространство. Плотности всех трех слоев для простоты положим одинаковыми. Скорость волны Рэлея для полупространства с параметрами каждого из трех слоев одинакова с точностью $ \sim {{10}^{{ - 5}}}$. Следовательно, дисперсия волны Рэлея в этом случае должна быть малой. Дисперсионные зависимости ${{V}_{R}}(\omega )$ и $E(\omega )$ показаны на рис. 9. Модель 3 представляет собой полупространство с параметрами первого слоя (табл. 1) и отвечает среде без дисперсии величин ${{V}_{R}}(\omega )$, $E(\omega )$.

Из представленных зависимостей следует, что относительное изменение ${{V}_{R}}(\omega )$ составляет около 1%, а относительное изменение $E(\omega )$ – около 3%. Таким образом, относительные изменения измеряемых величин действительно малы. При этом большие относительные изменения $E(\omega )$ указывают на информативность этого параметра. Отметим, что дисперсионные зависимости ${{V}_{R}}(\omega )$ для моделей 1 и 2 в полосе частот выше 10 Гц неразличимы, несмотря на существенное отличие коэффициентов Пуассона в этих моделях. Таким образом, имеется неоднозначность решения задачи определения параметров методами спектрального анализа поверхностных волн (SASW/MASW). В то же время частотная зависимость параметра эллиптичности $E(\omega )$ для моделей 1 и 2 отличается существенно в том же диапазоне частот. Следовательно, включение параметра $E(\omega )$ в анализ должно устранить неоднозначность определения параметров слоистой среды.

Сплошной тонкой линией показан результат численного моделирования (1) и (2) при уровне шума на 20 дБ ниже показанного на рис. 5, т.е. при использовании дополнительного когерентного накопления по ${{N}_{p}} = 100$ сеансам излучения. Результирующие погрешности, полученные при анализе данных в полосе 20–90 Гц, приведены в табл. 2. Нетрудно видеть, что, как и следовало ожидать, погрешности оценки максимальны для параметров полупространства. При этом погрешности оценки коэффициента Пуассона явно позволяют произвести диагностику насыщенности промежуточного слоя.

Таким образом, предложенное ранее [30, 32] расширение методов спектрального анализа поверхностных волн SASW/MASW за счет включения в анализ частотной зависимости параметра эллиптичности $E(\omega )$ действительно позволяет достичь точности, необходимой для решения важных задач дистанционной диагностики. Выражения (4) и (5), приведенные выше, носят общий характер и указывают путь достижения требуемой в каждом конкретном случае точности измерений. Хорошо известные выражения для когерентного временного (суммирование сеансов излучения) или пространственного (использование приемников в виде линейных антенн) накопления позволяют осуществить оценку необходимых ${{N}_{p}}$ и ${{N}_{g}}$.

Отметим еще одну интересную особенность совместного анализа ${{V}_{R}}(\omega )$ и $E(\omega )$. Выше отмечалось, что в стандартных хорошо известных приложениях SASW/MASW осуществляется анализ только дисперсионной зависимости ${{V}_{R}}(\omega )$. Это не позволяет получить оценки распределения скорости продольной волны и плотности по глубине. Добавление в анализ зависимости $E(\omega )$ снимает неоднозначность относительно скорости продольных волн. Плотность и ее зависимость от глубины в настоящей статье не рассматривается (параметр полагается известным). Здесь важно указать на существование достаточно точных эмпирических соотношений, связывающих плотность и скорость продольной волны в природных материалах, например, соотношения Гарднера $\rho \propto V_{P}^{{1/4}}$ и ему подобные [31]. Следовательно, имея возможность дистанционно оценивать скорости обеих объемных волн (${{V}_{P}}$ и ${{V}_{S}}$) на различных глубинах, мы получаем возможность дистанционной оценки плотности и ее зависимости от глубины. Точность такой оценки, очевидно, определяется не только точностью оценки ${{V}_{P}}$, но и погрешностями используемых эмпирических соотношений. Рассмотрение данного вопроса выходит за рамки настоящей статьи.

ВЫВОДЫ

Подведем итог выполненной работы.

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

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

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

4. Приведены примеры, имеющие тесную связь с ранее выполненными экспериментальными исследованиями. Это указывает на практическую реализуемость рассматриваемых в статье методов и дополнительно подтверждает полученные ранее экспериментальные результаты.

В сокращенном виде работа была представлена на XXXII сессии РАО в 2019 г. Исследования выполнены при поддержке гранта РФФИ 18-45-520006 и частичной поддержке госзадания ИПФ РАН по теме 0030-2021-0009.