Материал: Hydrogeodynamics101

До сих пор мы говорили об использовании интеграль­ного преобразования для получения аналитического или модельного решения той или иной краевой задачи. Одна­ко при исследовании некоторых вопросов фильтрации обратный переход от решения интегрального уравнения- аналога к решению исходного уравнения (от изображения к оригиналу) не является обязательным: искомые величи­ны могут определяться непосредственно из полученного решения в изображениях. Таковы, в частности, обратные задачи, связанные с определением фильтрационных пара­метров (см. гл. 5), в которых значения исходной функции S известны из полевых измерений. В этих случаях необ­ходимо предварительно рассчитать значения функции- изображения 15, используя известные из наблюдений гра­фики функции S(t). Для численного определения изобра­жения, отвечающего интегральному преобразованию Лапласа-Карсона вида (4.44), может использоваться сле­дующая приближенная формула [23]:

1 ⁰⁰ п

S = 7- / S{t) exp (- t/t\ dt = A₀ S(0) + 2 A, S(/_k).

^lp o ^v k —I

(4.53)

Коэффициенты А_к определяют по следующей табли­це.

Порядок вычислений изображения следующий: выби­рается значение параметра 1 (величина, имеющая раз­мерность времени), после чего из таблицы находят значе­ния и соответствующие им значения Далее, при известных значениях t и х_к, из соотношения t_K = х_к t_p определяют моменты времени t_K, на которые для вычис­лений по формуле (4.53) берут известные значения фун­кции S(t_K). Подсчет суммы ряда в формуле (4,53) обычно можно ограничить пятью-шестью первыми членами.

Если функция S(t) становится заметно отличной от нуля лишь при t > tj> 0, то вместо формулы (4.53) лучше использовать формулу

У = т* / S(t) e~^t/tP dt - е~^Р § A. S(t, + A t_k),

р о к — 1 (4.53 а)

где Л t_k = t_K — tj. В этом случае по таблице вместо t_K/t_p находим

Расчет изображения по графику функции S(f) занимает не­сколько минут.

Очевидно, что точность вычисления интегралов вида (4.44) в значительной степени связана с выбором параметра t . С одной сто­роны, величина t должна приниматься достаточно малой, т.е. зна­чение множителя exp(-t/t ) не должно быть слишком большим. Это положение определяется тем, что в выражении (4.44) интегрирова­ние по времени должно осуществляться в пределах от 0 до <», в то время как на практике фактические данные об изменении уровней подземных вод могут быть получены только в конечном интервале времени от 0 до t_Q.

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

В целом следует считать всегда желательным выполнение требо­вания:

t <— t

р ~ 6 ⁰ * <4.54а)

Если первые наблюдения до момента (_ш|п являются по каким-ли­бо причинам недостоверными, то следует принимать t отвечающим условию

*p^>2W (4.55)

Таким образом, для эффективного использования операционно­го метода должно выполняться условие (_ш5п 0,1 t_Q.

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

Рассмотрим, для примера, задачу об откачке из скважины в изолированном напорном пласте при произвольном дебите Q/0- Изображение для функции Q_c(t) определится формулой

Q_c=j-SQ_c(i)e-^,/,Pdt.

р о (4.56)

2ЛТ

Аналогично (4.50) решение поставленной задачи в изображени­ях дается формулой

(4.57)

илипри^г< 0,1*0,2

ЭД ₌_1 _1п 1,12

где

« Г ' (4.57а)

s(r) _ ^ко(хд

5_C ~K(xr_cy

Для совершенных скважин, работающих в условиях более слож­ных фильтрационных схем, решение (4.57) сохраняет свой вид, но коэффициент^ имеет отличные значения.

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

ПРИМЕР. Используем операционный метод для интерпретации режимных наблюдений, проведенных в паводковый период по створу пьезометров. Последние оборудованы на нижний слой в безнапорном двухслойном пласте. Створ ориентирован вкрест простирания реки, которая может считаться единственной гидродинамической грани­цей (полуограниченный пласт) и на которой задано условие третьего рода (3.56) (см. рис. 2.14).

Найдем сначала решение задачи в изображениях. Преобразуя исходное уравнение (4.6) по Лапласу-Карсону, получаем

3. Моделирование нестационарных плановых потоков

1. Конечно-разностная форма дифференциальных уравнений

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

с которым и пространство, и время разбиваются на конеч­ные интервалы, т.е. представляются дискретно на некото­рой пространственно-временной сетке с узловыми точка­ми x_j} у_;., t_k. При этом реальное непрерывное распределе­ние напоров Н(х, у, 0 заменяется дискретным: отыскива­ют или считают заданными напоры H(x_if y_t, t_k) во всех узловых точках сетки. Производные от напора в той или иной точке при этом заменяют приближенными конечно­разностными представлениями и они оказываются, таким образом, выраженными через разности в значениях напо­ров на концах пространственных или временных интерва­лов, включающих данную расчетную точку. Например, для одномерного плоскопараллельного случая, когда об­ласть фильтрации длиной L разбита сечениями x_t (0 < х,-< L) на интервалы длиной Дх, а расчетный период времени t разделен на последовательные промежутки t_A(0< t_k< t) продолжительностью Д t, имеем следующие выражения для производных в произвольной точке сетки Ц, /¹>:

д Н _ ^Я(*;, h) ~^н <Л-1> h) _ ^Я<Л+1> ⁽к)~^Н (^ХР h)

~ 9

дН

(4.64)

Дх

дх^{2 x}ih А*

Дх

Дх

(4.66)

(Д xf

В дальнейшем для упрощения записи мы введем след- щующую индексацию: H(x_it t_k) = н£ где н£— напор в

расчетном узле номер i (0< i<L/Ax на к-ом расчетном слое (0 £ к < A t/ t). В этих обозначениях, на­пример, уравнение (4.1) в точке /, к имеет следующее приближенное конечно-разностное представление:

_а* ^Яж - 2 Я,*+ яД_, ^ яД-Я*~'

“ (Дх)^{2 А<} ’ (4.67)

ЗАМЕЧАНИЕ. Если придать этому уравнению несколько иную форму:

яД-яД_, яД₊,-яД .яД-У-яД

Ах Ах Р At ’(4.67а)

то оно приобретает простой балансовый смысл (рис. 4.7,а). В самом деле, первое слагаемое в левой части - расход потока в среднем сечении (/ — 1/2) между узлами ini — 1, второе слагаемое — то же, в среднем сечении (/ + 1/2) между блоками / + 1 и /, а правая часть выра­жает собой скорость изменения объема воды, заключен­ного в интервале (/ — 1/2; i + 1/2), при снижении пьезо­метрической кривой. Иначе говоря, уравнение (4.67) — суть конечно-разностная форма условия сохранения мас­сы жидкости — уравнения неразрывности.

(4.60)

Общее решение этого обыкновенного дифференциального урав­нения с постоянными коэффициентами имеет вид [16 ]

А77 = Cj exp

Так как Л Н ^,=0, тоС₂ = 0. Значение С_} найдем из второго

граничного условия (при х^ш 0), которое в изображениях имеет вид (см. формулу (3.56)):

а (АН) _А#_р-А#_г

(4.62)

дх х—0 AL

где исходные функции-оригиналы АН_р и АН_г представляют собой

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

Д#= —7-г^Р— ехр

1+Д L^^^V\ Щ;)' (4.63)

Теперь, подсчитывая изображения от замеренных функций АН и АН , можно определить неизвестные параметры а и AL. Так,

(АЪ\

строя график связи In \~лг& для группы пьезометров, мы должны

получить прямую линию, угсш наклона СС которой к оси х дает значе­ние коэффициента уровнепроводности Vfl/ » ctg СС. Затем по отрез­ку Ь, отсекаемому на оси ординат, определяется параметр AL : 6 = ln(l + AL/\at\ Прямолинейность построенного гра­фика, которая должна наблюдаться для любых выбранных значений параметра преобразования является важным диагностическим признаком — свидетельством справедливости принятой расчетной схемы. При малом числе пьезометров (например, два) приходится ориентироваться на другой способ интерпретации. Сначала по отно-

АН_{ ^x2~~^xi шению = ехр ~ определяется коэффициент уровнепро­водности, а затем по известному значению АН для одного из пьезо­метров вычисляется параметр AJL Важный диагностический при­знак в этом случае — постоянство расчетных значений параметров при различных значениях параметра преобразования t.

2. Аналоговое моделирование нестационарной фильтрации

Придадим уравнению (4.67а) более общую форму,

считая, что параметры пласта ju* и Т могут меняться от сечения к сечению:

-t-k&t

y///////

LJULLL/JLLL

a

**4

a;

I.' ■| •

гггтгт

ax

71^г7~7~ГТТГ771~

L+i

6

Ui-tRi-i Ut *i« ^uh

Puc. 4.7. Моделирова­ние нестационарной фильтрации: a - исходная схема напорно­го пласта; 6 - схема рези­сторной сетки

где Т., и Т_м — проводимости на участках между сече­ниями i+(i - 1) и r*(rfl); вводя фильтра­ционные сопротивления согласно фор­муле (3.54)*, получаем:

н[L, -я? .Щн-н?

Ф,._ 1 ⁺ Ф,

* ¹

H_t (О.

'i+1 • • ^Ai ’ (4.68)

где w_i — площадь участка между сечениями i + 1 /2 и i - 1/2, прилежащего к точке i(a)_t. = Ах- 1);

И.

упругая водоотдача на этом участке.

Ширину потока В считаем равной единице.

В более компактной форме уравнение (4.68) можно записать, опуская индекс i расчетного узла:

v Hj—H * Н^к-Н^к~^х

Z ——ii ы гт—',

Ф} * A t (4.68а)

где индекс /=1,2 отвечает узлам, соседним с расчетным (я-2).

Рассмотрим теперь сетку электрических сопротивле­ний, элементы которой представлены на рис. 4.7б. Со­гласно закону Кирхгофа, сумма токов, поступающих в i-ый узел, должна равняться нулю:

U^k - U^k uL. - и^к 1/!~¹ - и^к

¹ _(_ ^{f 1 1} _|_ ‘ I _ Q

^Ri-1 ^Rt+1 ^R

t:

(4.69)

(обозначения потенциалов и и сопротивлений R ясны из рис. 4.76). Перепишем формулу (4.69):

Р*Ч -^и> vbi-ut_v}-dr

^Ri~ i ^Ri+1 ^RI, ’ (4.70)

или в более компактной форме:

и и^к- г/^-1

(4.70а)

Вспоминая теперь материал раздела 3.5.2, мы убежда­емся, что уравнения (4.70) и (4.70а) оказываются экви­валентными уравнениям (4.68) и (4.68а) нестационарной фильтрации, если потребовать соответствия:

a R - ^At

^ф ^ V®’ (4-71)

где а_ф — выбранный масштаб сопротивлений;

поэтому R_t носит название временного сопро­тивления.

Следовательно, на построенной таким образом сетке электрических сопротивлений потенциалы 17* отвечают

напорам Я* на расчетный момент &, — если на концы временных сопротивлений R_t подаются известные потен-