Это предварительная версия научной статьи. После публикации здесь будут размещены DOI, выходные данные и корректная ссылка для цитирования.
Не используйте этот препринт для научного цитирования.
Рассмотрена задача перехода от литологических колонок буровых скважин к трёхмерной цифровой инженерно-геологической модели. Выполнен обзор существующих подходов к выделению инженерно-геологических элементов и к восстановлению их пространственного распространения: нормативной схематизации по коэффициенту вариации, кластерного и дискриминантного анализа, классификации по параметрам статического зондирования, явных и неявных методов построения границ, категориальных геостатистических методов и марковских моделей переходных вероятностей. Показано, что явные методы, интерполирующие отметки кровель и подошв, принципиально не воспроизводят переслаивание, при котором один элемент входит в колонку многократно. Предложена методика, в которой интерполируется не характеристика и не отметка границы, а принадлежность точки пространства элементу: строятся индикаторы элементов, выполняется ординарный кригинг индикаторов с общей матрицей системы, вероятности нормируются, а мерой локальной неопределённости служит информационная энтропия. Методика реализована в программе KrigingBH и проверена на площадке из шестнадцати скважин с пятью элементами, в четырнадцати из которых зафиксировано переслаивание. Индикаторный метод воспроизвёл литологические колонки в собственных узлах на 97,6 % против 87,4 % у явного метода; при перекрёстной проверке с исключением выработок совпадение составило 70,4 % при среднеквадратической ошибке отметок границ 1,38 м. Показано, что два приёмочных критерия — доля объёма с высокой энтропией и точность положения границ — дают противоположные заключения, и предложено контролировать оба.
Переход от бумажных разрезов к цифровой инженерно-геологической модели (ЦИГМ) поставил задачу, которая при ручном построении решалась инженерным суждением и потому не формализовалась: как по конечному числу колонок восстановить пространственное распространение инженерно-геологических элементов (ИГЭ) и как оценить достоверность полученного результата.
Задача распадается на две части. Первая — схематизация: объединение разностей грунтов в инженерно-геологические элементы (далее в тексте — элементы), для которых допустимо принимать единые нормативные и расчётные значения характеристик. Вторая — оконтуривание: построение границ выделенных элементов между выработками. В нормативных документах формализована лишь первая часть, причём через статистический критерий однородности; вторая оставлена на усмотрение исполнителя. Между тем именно она определяет вид разреза, объёмы элементов и, в конечном счёте, расчётные результаты.
Настоящая работа посвящена методике, в которой обе части увязаны, а результат сопровождается количественной оценкой неопределённости. Изложению методики предпослан обзор существующих подходов, необходимый для того, чтобы показать, какие именно затруднения она устраняет.
Нормативная схематизация опирается на проверку однородности выборки по коэффициенту вариации:
где X — характеристика грунта; n — число единичных определений; [V] — допускаемое значение, различное для разных характеристик [1, 2].
Элемент считается выделенным правильно, если для всех контролируемых характеристик условие выполняется; в противном случае элемент подлежит разделению. Подход прост, встроен в практику и обеспечивает связь с расчётными значениями, однако имеет два ограничения. Во-первых, он оперирует выборкой, уже отнесённой к элементу, то есть проверяет гипотезу, но не порождает её: само отнесение выполняется по геологическим соображениям. Во-вторых, он не содержит пространственной составляющей — две пробы из противоположных концов площадки входят в одну выборку наравне с соседними.
Попытки формализовать подобный подход привели к применению методов многомерной классификации. Кластерный анализ по совокупности характеристик позволяет выделять группы без предварительного задания их числа [3, 4]; дискриминантный анализ решает обратную задачу отнесения новой пробы к уже выделенным группам [5]. Развитием служат методы машинного обучения — от деревьев решений до искусственных нейронных сетей, применяемые в основном к данным непрерывного зондирования [6, 7]. Общая слабость этих методов в приложении к изысканиям состоит в том, что они игнорируют стратиграфическую и генетическую принадлежность: статистически близкие, но разновозрастные грунты попадают в один класс, что противоречит требованиям нормативных документов [2].
Отдельную линию образует классификация по параметрам статического зондирования через индекс типа поведения грунта [8, 9]. Она даёт непрерывное по глубине разделение и потому особенно удобна для построения моделей, но, как неоднократно отмечалось [10], не является определением вида грунта по классификации [11] и требует тарирования по совмещённым скважинам.
Методы построения границ удобно разделить на три группы.
Явные методы интерполируют отметки кровель и подошв элементов как непрерывные поверхности:
где Zk — отметка границы k-го элемента; λ — веса интерполятора, в качестве которого применяют кригинг, метод обратных расстояний, сплайны или радиальные базисные функции [12, 13].
Подход даёт гладкие поверхности, удобен для передачи в системы проектирования и прозрачен для проверки. Его принципиальное ограничение состоит в предположении, что каждый элемент входит в колонку однократно и что порядок элементов сверху вниз одинаков во всех выработках. Переслаивание, линзы и выклинивание в такую схему не укладываются: при многократном вхождении элемента приходится либо объединять прослои, либо вводить искусственные поверхности, число которых заранее неизвестно.
Неявные методы описывают границу как множество уровня скалярного поля, построенного по данным о принадлежности и по ориентировкам контактов [14, 15]. Они снимают ограничение на топологию тел и широко применяются в рудной геологии. Затруднение при переносе в изыскания связано с исходными данными: поле строится по контактам с известной ориентировкой, тогда как вертикальная скважина даёт только отметку пересечения, а падение слоя остаётся неизвестным.
Категориальные геостатистические методы работают непосредственно с признаком принадлежности. Индикаторный кригинг [16, 17] даёт оценку вероятности принадлежности каждому классу; последовательное индикаторное моделирование [18] порождает ансамбль равновероятных реализаций, сохраняющих заданную структуру; статистика многоточечных шаблонов [19] воспроизводит формы тел, заданные обучающим образом. Отдельную группу составляют марковские модели: цепи Маркова в вертикальном направлении и модели переходных вероятностей [20, 21], в которых структурная информация задаётся матрицей переходов, имеющей прямой геологический смысл — средние мощности и последовательность напластования. Байесовские постановки [22, 23] объединяют оценку числа слоёв, положения границ и характеристик в единой вероятностной схеме.
Сопоставление методов по признакам, существенным для изыскательской практики, приведено в табл. 1.
Таблица 1. Сопоставление подходов к восстановлению распространения ИГЭ
Из обзора следует, что для изыскательской практики предпочтительны категориальные методы, а среди них — индикаторный кригинг как наименее требовательный к исходным данным и вычислительно наиболее дешёвый. Однако его применение в изысканиях сдерживают три обстоятельства: необходимость строить и обосновывать вариограмму для каждого элемента отдельно; нарушение порядковых соотношений, вследствие которого оценки вероятностей выходят за пределы отрезка от нуля до единицы и не суммируются в единицу; отсутствие принятого способа предъявить полученную неопределённость в форме, пригодной для приёмки модели.
Целью работы является методика, устраняющая эти три проблемы и доведённая до проверяемого расчётного результата.
Последовательность операций приведена на рис. 1. Исходными данными служат литологические колонки выработок с выделенными по (1) элементами и таблица расчётных значений характеристик по элементам. Существенно, что схематизация выполняется до построения модели и средствами методики не изменяется: модель восстанавливает распространение уже выделенных элементов, а не подменяет инженерное суждение.
Рис. 1. Последовательность операций методики
Колонки дискретизируются с шагом, не превышающим половины наименьшей мощности прослоя, подлежащего воспроизведению. Для каждого элемента строится индикатор (3). Структурный анализ выполняется по объединённой индикаторной вариограмме (7), после чего единожды собирается и разлагается матрица системы (4), а решения находятся по числу элементов. Вероятности нормируются по (8), наиболее вероятный элемент определяется по (9), неопределённость характеризуется энтропией (10). Верификация выполняется перекрёстной проверкой с исключением выработок и оценкой доли объёма с высокой энтропией.
Принадлежность точки элементу описывается индикатором:
Оценка выполняется ординарным кригингом. Система уравнений в матричной форме
где Γ — матрица значений вариограммы между опорными точками; λ — параметр регуляризации; E — единичная матрица; ω и μ — веса и множитель Лагранжа; i — вектор значений индикатора.
В данной работе применена сферическая модель вариограммы:
с остаточной дисперсией С0, порогом C и радиусом a.
Для индикаторов порог принимается равным единице: веса ординарного кригинга инвариантны к масштабу порога, поэтому его абсолютное значение на результат не влияет, а остаточная дисперсия и регуляризация задаются долями порога. Геометрическая анизотропия учитывается приведением лага:
Первая из трёх проблем, отмеченная в п. 2.3, устраняется структурой системы: матрица (4) зависит только от расположения опорных точек и модели вариограммы, но не от номера элемента. Следовательно, разложение выполняется однократно, а решения по элементам получаются обратным ходом. Вычислительные затраты составляют порядка n³/3 операций на разложение и K·n² на решения вместо K·n³/3 при поэлементном построении, что при пяти элементах даёт сокращение почти впятеро.
Обоснование модели по вариограммам отдельных индикаторов требует K независимых подгонок, причём для малопредставленных элементов оценка оказывается неустойчивой. Предлагается использовать объединённую индикаторную вариограмму:
равную вероятности того, что две точки, разделённые вектором h, принадлежат разным элементам. Величина вычисляется по всем парам опорных точек одновременно, имеет прямой геологический смысл — вероятность смены элемента на заданном расстоянии — и изменяется от нуля при нулевом лаге до значения, определяемого долями элементов, при больших лагах. Направленные оценки по вертикали и по горизонтали дают отношение радиусов, принимаемое за отношение анизотропии.
Индикаторный кригинг не гарантирует принадлежности оценок отрезку от нуля до единицы и их суммирования в единицу. Второе затруднение устраняется усечением и нормировкой:
после чего строение определяется правилом максимума:
Третье затруднение — предъявление неопределённости — решается двумя показателями. Локальная неопределённость характеризуется информационной энтропией и её нормированным значением
Энтропия равна нулю, если элемент определён однозначно, и достигает ln K при равновероятности; выражена в натах, поскольку используется натуральный логарифм. Интегральной характеристикой служит доля объёма модели, в которой нормированная энтропия превышает заданный порог:
Эта величина пригодна в качестве приёмочного критерия: она не зависит от числа элементов и допускает нормирование.
Достоверность оценивается перекрёстной проверкой с последовательным исключением выработок: выработка удаляется из набора, модель перестраивается по остальным, и её колонка в точке исключённой выработки сопоставляется с фактической. Вычисляются доля совпадения элементов и оценка Брайера:
а также среднеквадратическая ошибка положения границ, определяемая как расстояние от каждой фактической границы до ближайшей модельной границы:
Первые два показателя характеризуют качество вероятностного прогноза, третий — точность геометрии. Как показано далее, они не эквивалентны и могут давать противоположные заключения.
Расчётные значения характеристик присваиваются выделенным элементам, после чего колонка модели в произвольной точке используется в расчёте осадки методом послойного суммирования [24]:
где модуль деформации берётся по наиболее вероятному элементу в середине элементарного слоя. Тем самым расчёт выполняется не по осреднённому разрезу, а по колонке модели в точке фундамента.
Методика реализована в программе KrigingBH [25] и применена к площадке, представленной шестнадцатью буровыми скважинами (рис. 2). Общий объём данных составляет 96 слоёв, пять инженерно-геологических элементов, три стратиграфо-генетических подразделения. Площадка имеет размеры 213 × 77 м, глубина скважин 15,0 и 25,0 м, расстояние между соседними выработками от 22,3 до 56,0 м при среднем 27,4 м. Категория сложности инженерно-геологических условий вторая.
Существенная особенность набора: в четырнадцати выработках из шестнадцати зафиксировано переслаивание, при котором один элемент входит в колонку многократно — до трёх раз. Мощности отдельных прослоев составляют 0,8—2,2 м. Именно этот набор делает сопоставление подходов показательным.
Рис. 2. Ситуационный план: шестнадцать буровых скважин и линия продольного разреза через восемь выработок южного ряда
Дискретизация с шагом 0,5 м дала 760 опорных точек. Объединённые индикаторные вариограммы приведены на рис. 3. В горизонтальном направлении вариограмма растёт медленно, от 0,29 на лаге 27 м до 0,40 на лаге 135 м; формальная подгонка даёт радиус 239 м, превышающий размеры площадки, то есть горизонтальный радиус корреляции по имеющейся сети не определяется, а лишь ограничивается снизу. В вертикальном направлении вариограмма уверенно выходит на плато: от 0,16 на лаге 0,7 м до 0,93 на лаге 12,5 м при подобранном радиусе 13,0 м.
Для построения приняты радиус 150 м и отношение анизотропии 20, что отвечает вертикальному радиусу 7,5 м, остаточная дисперсия 0,10 от порога и параметр регуляризации 0,01. Вертикальный радиус принят вдвое меньше подобранного сознательно: подгонка выполнена по данным, включающим вертикальный тренд слоистости, и завышает радиус, тогда как мощности переслаивающихся прослоев составляют порядка метра.
Рис. 3. Объединённая индикаторная вариограмма по 760 опорным точкам: горизонтальное (а) и вертикальное (б) направления
Продольный разрез построен по ломаной линии, длиной 227 м через восемь выработок на сетке 280 × 160 узлов (рис. 4). Модель воспроизводит трёхчленное строение: перекрывающий насыпной грунт мощностью 1,0-3,4 м, толщу переслаивания песка и суглинка в интервале отметок 135-146 м и подстилающие глины ниже отметки 138 м.
Поле энтропии локализует неопределённость: она сосредоточена в узких приграничных зонах и в интервале переслаивания, где достигает 1,2-1,5 нат при максимуме ln 5 = 1,609 нат, тогда как в толще нижних глин не превышает 0,3 нат. Средняя по объёму энтропия составила 0,447 нат, доля объёма с нормированной энтропией выше 0,8 — 0,64 %.
Рис. 4. Продольный разрез: поле наиболее вероятного элемента с наложением фактических колонок (а) и поле информационной энтропии (б)
Объёмное поле построено на сетке 56 × 40 × 56, что составляет 125 440 узлов; расчёт занимает около секунды. Доли объёма по элементам: 11,0; 18,3; 11,2; 7,3 и 52,1 %. Разбиение поля на связные тела дало семь тел при этом, третий и четвёртый элементы распались на два тела каждый, что отвечает переслаиванию, — и 43 332 треугольника замкнутых оболочек, передаваемых в обменный файл IFC 4.3 вместе с показателями неопределённости (рис. 5).
Рис. 5. Трёхмерная модель: объём элементов с вырезом ближнего угла (а) и поле модуля деформации как математическое ожидание по вероятностям (б)
Перекрёстная проверка по шестнадцати выработкам занимает около пяти секунд. Результаты сведены в таблице 2.
Таблица 2. Результаты верификации
Разброс по выработкам закономерен. Наилучшие показатели получены для выработок, расположенных внутри сети и вскрывающих простое трёхчленное строение; наихудший результат, у краевой выработки, колонка которой восстанавливается экстраполяцией и в которой суглинок залегает непосредственно под насыпным грунтом мощным слоем 7,4 м, тогда как в соседних выработках эту позицию занимает песок.
Для условного фундамента размером в плане 3 × 3 м с глубиной заложения 4,0 м при среднем давлении 300 кПа выполнен расчёт осадки по колонке модели в точке каждой из шестнадцати выработок (рис. 6). Осадка составила от 44,2 до 56,3 мм при среднем 51,5 мм и коэффициенте вариации 0,073, мощность сжимаемой толщи 4,0-4,5 м.
Малый разброс при среднеквадратической ошибке положения границ 1,38 м объясняется тем, что сжимаемая толща при указанной геометрии целиком располагается в интервале переслаивания, модули деформации элементов которого различаются лишь на четверть; перестановка прослоев внутри толщи почти не меняет суммы (14). Вывод не переносится на другие задачи: для свайных фундаментов, где определяющим является положение кровли несущего слоя, та же ошибка дала бы существенно больший разброс.
Рис. 6. Сопоставление фактической колонки скважины с колонкой модели в той же точке (а) и расчётная осадка по колонкам модели во всех шестнадцати выработках (б)
Для количественного сопоставления по тем же данным построен разрез явным методом: для каждого элемента взята отметка кровли первого вхождения, отметки интерполированы в плане ординарным кригингом по (2) с последующим приведением к монотонной последовательности сверху вниз. Результат приведён на рис. 7 в сопоставлении с индикаторным построением по той же линии.
Рис. 7. Продольный разрез, построенный явным методом (а) и индикаторным кригингом (б); столбиками показаны фактические колонки выработок
Явный метод даёт слоистую картину, в которой переслаивание песка и суглинка объединено в два монолитных слоя, а четвёртый элемент отсутствует: его кровля при приведении к монотонной последовательности оказывается ниже кровли пятого элемента. Совпадение с фактическими колонками в собственных узлах выработок составило 87,4 % против 97,6 % у индикаторного метода, причём в краевой выработке с мощным суглинком оно уменьшилось до 38 %.
Следует от метить, что различие носит не количественный, а качественный характер. Явный метод не способен воспроизвести многократное вхождение элемента в литологическую колонку, так как число интерполируемых поверхностей задаётся заранее и равно числу элементов. Индикаторная постановка такого ограничения не содержит, поскольку интерполируется признак принадлежности в точке, а число и форма тел получаются как следствие.
Наиболее существенным для практической реализации представляется результат сопоставления двух приёмочных критериев. Критерий по энтропии выполняется с большим запасом: доля объёма с высокой неопределённостью составляет 0,64% при типичном допуске 15%. Критерий по среднеквадратической ошибке отметок границ, напротив, не выполняется: полученные 1,38 м превышают допуск, разумный для второй категории сложности инженерно-геологических условий.
Противоречие объясняется тем, что показатели измеряют разное значение. Энтропия оценивает долю объёма, в которой модель не уверена, и эта доля мала, поскольку высокая энтропия сосредоточена в тонких приграничных зонах. Ошибка отметок границ оценивает точность положения самих границ, и она велика, поскольку мощности переслаивающихся прослоев сопоставимы с произведением расстояния между выработками на отношение анизотропии. Отсюда следует, что энтропийный критерий необходим, но недостаточен, и в состав обязательных метрик приёмки модели должны входить обе величины.
Второе замечание касается предела достижимой точности. При шаге сети 22-56 м и мощностях прослоев порядка метра положение границ переслаивающихся тел между выработками достоверно не определяется никаким методом интерполяции: информации для этого в данных нет. Значение методики в этом случае состоит не в том, чтобы дать точные границы, а в том, чтобы честно показать, где их положение неизвестно, и не создавать ложного впечатления определённости, которое возникает при построении гладких поверхностей явным методом.
Третье замечание относится к схематизации. Методика восстанавливает распространение уже выделенных элементов и не подменяет их выделение. При этом результат верификации может служить обратной связью: устойчиво низкое совпадение по отдельным выработкам указывает либо на недостаточность сети, либо на то, что схематизация объединила разнородные тела.
Предложена методика построения цифровой инженерно-геологической модели, в которой интерполируется принадлежность точки пространства инженерно-геологическому элементу. Три проблемы, сдерживавшие применение индикаторного кригинга в изысканиях, устранены: общая для всех элементов матрица системы сокращает вычислительные затраты почти впятеро при пяти элементах; объединённая индикаторная вариограмма заменяет множество независимых подгонок одной, имеющей прямой геологический смысл; усечение с нормировкой восстанавливает порядковые соотношения вероятностей, а энтропия и доля объёма с высокой энтропией дают показатели неопределённости, пригодные для приёмки.
На площадке из шестнадцати скважин с переслаиванием в четырнадцати из них методика воспроизвела литологические колонки в собственных узлах на 97,6 % против 87,4 % у явного метода интерполяции отметок границ, причём явный метод не воспроизвёл переслаивание и потерял один из пяти элементов. При перекрёстной проверке с исключением выработок совпадение составило 70,4 % при среднеквадратической ошибке отметок границ 1,38 м.
Показано, что два приёмочных критерия, такие как доля объёма с высокой энтропией и точность положения границ, могут давать противоположные заключения, и в этом случае предлагается контролировать их совместно. Показано также, что при сжимаемой толще, полностью расположенной в интервале переслаивания близких по деформируемости элементов, пространственная изменчивость расчётной осадки фундамента незначительна (коэффициент вариации 0,073) несмотря на значительную ошибку положения границ.
Дальнейшее развитие методики видится в совместном использовании данных бурения и статического зондирования: калибровка по совмещённым парам позволяет привлечь непрерывные записи зондирования как дополнительные опорные точки, что должно существенно улучшить разрешение по вертикали при сохранении нормативной схематизации по данным бурения.
наверх
Поделиться
Поделиться
Скопировать ссылку
Скопировать ссылку 
Vk 
Telegram 
Viber 
WhatsApp 
Оформить подписку