Статистическое описание осциллятора с учетом флуктуирующего коэффициента трения Текст научной статьи по специальности «Физика»

ФИЗИКА

J

УДК 519.62

А. Н. Морозов, А. В. Скрипкин

СТАТИСТИЧЕСКОЕ ОПИСАНИЕ ОСЦИЛЛЯТОРА С УЧЕТОМ ФЛУКТУИРУЮЩЕГО КОЭФФИЦИЕНТА ТРЕНИЯ

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

Традиционное описание необратимых процессов, в частности броуновского движения, не учитывает экспериментально наблюдаемых флуктуаций кинетических коэффициентов, которые имеют спектр типа фликкер-шума [1, 2]. Для описания процессов со спектром этого типа неприменимы подходы, основанные на теории марковских случайных процессов [3-5], так как с помощью конечной системы дифференциальных уравнений невозможно преобразовать процесс с независимыми приращениями в процесс типа фликкер-шум.

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

Уравнение Ланжевена. Уравнение движения броуновской частицы в вязкой среде можно записать в виде

Р + з = г (г) + $ (г), (1)

где Р — импульс броуновской частицы; 3 — термодинамический поток; Г (г) — внешняя детерминированная сила; $ (г) — случайное воздействие. В линейном приближении термодинамики необратимых процессов поток 3 связан с термодинамической силой X следующим соотношением [6]:

г

з (г) = ! с (г, т) х (т) ¿т, (2)

0

где С (£, т) — ядро преобразования, которое обычно представляется в виде

С (£,т) = 2Ш (£ - т). (3)

Здесь Д — коэффициент диффузии, который считается постоянным.

Термодинамическая сила X для броуновской частицы определяется по формуле

х = м£т, (4)

где М — масса броуновской частицы; — постоянная Больцмана; Т — температура среды.

Если ядро преобразования С (£, т) в форме (3) подставить в формулу (2), а затем в уравнение (1), то оно приобретает вид уравнения Ланжевена [7, 8], которое используют для описания броуновского движения:

Р + пР = F (£) + £ (£). (5)

Здесь п — коэффициент вязкого трения, определяемый по формуле

п = мквт. (6)

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

Отметим также, что процесс Р (£), описываемый уравнением (5), представляет собой марковский случайный процесс, если случайное воздействие £ (£) есть среднеквадратическая производная процесса с независимыми приращениями, например винеровского процесса.

Уравнение броуновского движения, учитывающее флуктуации коэффициента трения. Для описания броуновского движения с учетом флуктуаций коэффициента вязкого трения предположим, что ядро преобразования С (£, т) является функцией, изменяющейся случайным образом в зависимости от значений £ и т. Ядро преобразования С (£, т) в общем случае представимо в виде [6]

С (£,т) = (£1 (£ - т) £ (т)>, (7)

где (.. .> — операция нахождения математического ожидания; £ (т) — среднеквадратическая производная процесса с независимыми приращениями W (т); £1 (£ — т) — средняя случайная сила на интервале

времени t — т, т.е.

6 (t — т) = . (8)

t — т

Дальнейший интерес представляет случай, когда G (t, т) случайным образом зависит от значений t и т, поэтому без учета оператора нахождения математического ожидания вместо (7) имеем

G (,т ) = W (t — т) d^ir), (9)

t — т ат

где W (t) — процесс с независимыми приращениями, который связан со случайным воздействием £ (t) интегралом Ито:

t

w (t) = / £ (т) ат. (10)

0

В этом случае выражение (2), связывающее поток J с силой X, принимает вид

t

J (t) = j w,^X (т) ат. (11)

0

Подстановка этой формулы в выражение (1) позволяет вместо уравнения Ланжевена (5) получить следующее уравнение для описания броуновского движения:

t

• [W (t — т) dW (т) P (т) ,

P +j ^—-Т-ОТ2Mpk(iОт = F (t) + £ (t), (12)

0

учитывающее флуктуации коэффициента трения броуновской частицы.

Случайный процесс P (t), описываемый стохастическим интегральным уравнением (12), является немарковским процессом, и для определения его статистических характеристик неприменима теория стохастических дифференциальных систем [8]. Из-за случайного характера изменения импульса P (т) для такого уравнения неприменима и теория линейных немарковских процессов [9] и требуется численное моделирование.

Спектральная плотность флуктуаций кинетических коэффициентов. Рассмотрим случай, когда термодинамическая сила X представляет собой постоянную величину: X = const. Тогда соотношение (2) может быть записано в виде

J (t) = y (t)X, (13)

где кинетический коэффициент 7 рассматривается как случайный процесс, описываемый выражением

г Г

7 (()_ { ^ = { ^^ (т). (14)

0 0

Далее полагаем, что процесс И (г) представляет собой винеров-ский случайный процесс с интенсивностью V = 2Д, связанный со случайным воздействием £ (г) формулой (10). Тогда математическое ожидание кинетического коэффициента 7 (г) равно

t . t (7(t)> (t - T > dW (t )\ = /■ (W (t - T >(T»

г - т 4 ' / у г - т

г , 2\ г 0 (15)

Г (.И (т))2 Г ¿т г

= --V -« Vш—,

У г - т У г - Т гг'

00

где гг — малая положительная величина. При нахождении интеграла в этом выражении учтено, что процесс И (г) является процессом с независимыми приращениями и, следовательно, для него справедливо соотношение

((¿И )2) = ^г. (16)

Из выражения (15) следует, что при описании кинетических коэффициентов в рамках рассматриваемой модели их средние значения возрастают с течением времени по логарифмическому закону. Малый промежуток времени гг можно считать величиной, близкой к значению постоянной времени хаотизации частиц рассматриваемой системы. В частности, для разреженного газа величина гг равна среднему времени соударения его молекул [7].

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

Определим момент второго порядка для процесса 7 (г):

(7(«1)7(«2» = 1 ¿И (г,)) ¿И Ы> =

г1 - Т1 / г2 - т2

t2 tX (W (ti - Ti) dW (ti) W (t2 - T2) dW (T2)>

(¿1 - Ti) (¿2 - T2) 0 0

t2 1i

2 [ f (1 + 2(min (¿1 ,t2) - max(Ti,T2) 8 (т1 - T2))) , ,

= ' J J -(tl-Ti)(t2 -T2)-'Т2'Т1 =

0 0

t2 ti

2 (1 + 2 (¿1 - max (ti, T2) 8 (ti - T2)))

= v II -т.-гтг-ч-dT2dT1 =

(¿1 - T1) (¿2 - T2)

00

= Л 1П £ 1П + 2 1П-^ ), (17)

V & ¿2 - ¿1 + &)

где принято ¿2 ^ ¿1.

Из формулы (17) следует выражение для ковариационной функции

¿2

K (¿1, ¿2) = (7 (¿1) 7 (¿2)) - (7 (¿1)) (7 (¿2)) = 2V ln , , + 8,

t2 - t1 + 8t (18)

которая в свою очередь позволяет определить одностороннюю спектральную плотность случайного процесса 7 (¿):

( 7' \

G (ш) = lim lim < 4 K (£ - т, ¿) cos шт^т > =

х 7 t^^ ät^o 1 J

= lim < 8v 1 ln — cos шт^т > =-. (19)

Ч-го I J Tш

Таким образом, флуктуации кинетического коэффициента, описываемого выражением (14), имеют спектральную плотность

4nv

G (ш) = —, (20)

ш

обратно пропорционально зависящую от частоты. Это указывает на то, что флуктуации кинетического коэффициента 7 (¿) представляют собой фликкер-шум.

Математическое моделирование движения броуновской частицы под действием постоянной силы. Численное решение уравнения (12) позволило установить особенности движения броуновской частицы в среде с флуктуирующим коэффициентом вязкого трения. При численном решении этого уравнения считалось, что интенсивность винеровского процесса W (¿) равна v = 2 (дисперсия D = 1), произведение MkBT = 1, а шаг по времени Ai = 0,3. Значение шага Ai определено из условия устойчивости численного решения. Внешнюю силу варьировали в диапазоне F = 0,00... 1,00 с шагом AF = 0,05. Всего в расчетах для одного значения силы было получено 105 зна-

Рис. 1. Зависимость математического ожидания Мр (1) и дисперсии Ор (2) от внешней силы Р

чений импульса Р (£), по которым определялись математическое ожидание Мр, дисперсия Бр и спектральная плотность Ср (ш) процесса

р (*).

На рис. 1 приведены зависимости математического ожидания Мр и дисперсии Вр флуктуаций импульса броуновской частицы Р (£) от внешней приложенной силы Р. Видно, что зависимость Мр (Р) имеет линейный характер с коэффициентом пропорциональности, близким к единице, что совпадает с традиционным описанием при п =1. Зависимость дисперсии флуктуаций импульса Рр (Р) близка к квадратичной, что связано с приближением средней скорости движения броуновской частицы к средней квадратичной скорости ее хаотического движения.

На рис. 2 приведена спектральная плотность Ср (ш) флуктуаций импульса Р (£) для единичной внешней силы (Р = 1). Спектр рассчитан по 2-106 точкам. Из полученной зависимости Ср (ш) следует, что

Gp( ш)

0,0008

0,0007 0,0006 0,0005 0,0004 0,0003 0,0002 0,0001 0

1

0,0002 0,0004 0,0006 0,0008 Ю

Рис. 2. Спектральная плотность флуктуаций импульса Р(£) при внешней силе Р = 1

при действии внешней силы флуктуации импульса броуновской частицы приобретают характер фликкер-шума. Аппроксимация для первых двухсот точек дает зависимость Ор (ш), близкую к гиперболической:

СРН = 7 • 10-9ш-0'912. (21)

Таким образом, описание броуновского движения в среде с флуктуирующим коэффициентом вязкого трения как немарковского случайного процесса более реалистично, чем традиционное, с использованием уравнения Ланжевена. При этом флуктуации импульса броуновской частицы содержат не только составляющие, описываемые белым шумом, но и флуктуации со спектром типа фликкер-шума, возникающие при приложении к частице внешнего детерминированного воздействия Р (¿).

Броуновская частица как осциллятор. Рассмотрим случай, когда на броуновскую частицу помимо детерминированной силы Р (¿) и случайного воздействия £ (¿) действует возвращающая сила, пропорциональная отклонению частицы У (¿) от положения равновесия,

^ (*) = —кУ (*), (22)

где к — коэффициент, характеризующий упругие свойства среды, не зависящий от времени. В этом случае уравнение Ланжевена (5) имеет вид

Р + пР + кУ (¿) = Р (¿) + £ (¿). (23)

Из уравнения (23) легко может быть найдена спектральная плотность флуктуаций импульса частицы Р(¿). Для этого применим преобразование Лапласа к выражению (23), считая силу Р(¿) равной нулю, а случайную силу £(*) — белым шумом. Получим

вР(в) + пР(з) + кМ^ = ¿(Я), (24)

где в — параметр преобразования, Р(в) — изображение функции Р(¿),

£(в) — изображение белого шума. Из выражения (24) находим

в

Р(в) = 2 + + 2 £(в), (25)

в2 + пв + ш0

где ш0 = \J~~k — собственная частота осциллятора. Учитывая, что спектральная плотность

2

(26)

а спектральная плотность белого шума, равная его интенсивности, —

постоянная величина (а = D), окончательно получим

w 2

°Р(") = 4 \ ( 2 О 2 \ 4 V- (27)

и4 + (п2 - 2^2)ч2 + и)

Броуновская частица как осциллятор в среде с флуктуирующим коэффициентом трения. Если осциллятор находится в среде с флуктуирующим коэффициентом трения, то вместо дифференциального уравнения Ланжевена (23) для его описания применимо уравнение

г

р + 1 ^ ^^ = м + р т ((,). (28)

0

Уравнения (28) решается численно. При вычислениях шаг по времени принят Д£ = 0,3, произведение МквТ = 1, М = 1, интенсивность белого шума V = 1. Внешняя сила варьировалась в диапазоне Р = 0,0 ... 1,0 с шагом ДР = 0,1. Параметр к принят равным 0,1. Кроме того, расчеты проводили при заданных значениях внешней силы (Р = 0,0 и Р = 0,5) и разных значениях к в диапазоне к = 0,05 ... 0,3 с шагом Дк = 0,05. Все величины формально считали безразмерными. В каждом случае получено 105 значений импульса и координаты осциллятора.

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

Му = 9,932Р - 0,081. (29)

1 / ---- у 2

/у ' у / / / X /

у/ X / / У

У У

/X у/ у / у /

0 0,2 0,4 0,6 0,8 ^

Рис.3. Зависимость математического ожидания координаты частицы Му от приложенной внешней силы Р (1) и прямая аппроксимации (2)

40 35 30 25 20 15 10

Dy(k)

у

\ \

\ \ \ \ \ \

v V \\

2 ч\

\

NN

0,05

ОД

0,15

0,2

0,25

При увеличении силы на Д^ = 0,1 среднее значение координаты частицы с большой степенью точности увеличивается на ДУ = 1. Такое отклонение связано с выбором к, так как из классического описания осциллятора следует, что ДУ = Д^/к. Полученное совпадение с результатом классической теории служит одним из подтверждений адекватности используемой модели реальному процессу.

На рис.4 изображена зависимость дисперсии флуктуаций координаты Бу от коэффициента упругости к при ^ = 0. Видно, что дисперсия уменьшается нелинейно с ростом к. Здесь же приведена аппроксимирующая функция, соответствующая степенному закону (штриховая линия). Видно, что дисперсия флуктуаций координаты уменьшается обратно пропорционально коэффициенту к. Полученное с помощью метода наименьших квадратов уравнение аппроксимирующей кривой имеет вид

Бу = 0,934к-1'192. (30)

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

Рис. 4. Зависимость дисперсии флуктуаций координаты от коэффициента упругости к при Р = 0 (1) и степенная функция аппроксимации (2)

1 k (У2> = 1 kBT.

(31)

2 х ' 2

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

Бу = <У2) . (32)

Выражение (32) следует из определения дисперсии Бу = (У2) — (У )2 с учетом того, что для осциллятора в отсутствие внешней силы (У) = 0. Из соотношений (31) и (32) следует

кв Т

Dy =

k

(33)

Таким образом, зависимость дисперсии координаты классического осциллятора Бу от коэффициента упругости к соответствует обратно

пропорциональной зависимости, что подтверждается результатами численного расчета (см. формулу (30)).

Дисперсия флуктуаций импульса осциллятора Вр почти не меняется с ростом коэффициента упругости к (соответствующая зависимость приведена на рис. 5) и принимает значение, близкое к единице. Это также согласуется с классической моделью. Действительно, так как

Рис. 5. Зависимость дисперсии флуктуаций импульса Ор от коэффициента упругости к при Р = 0

2M 2> = 2т,

а дисперсия импульса определяется выражением

DP = (P2> ,

(34)

(35)

находим

Бр = Мкв Т, (36)

т.е. дисперсия флуктуаций импульса зависит только от массы осциллятора и температуры среды и не зависит от коэффициента упругости. Незначительный рост дисперсии флуктуаций импульса осциллятора Вр с увеличением параметра упругости к (см. рис. 5) связан, по-видимому, с возрастанием ошибок численного моделирования при увеличении собственной частоты колебаний осциллятора.

На рис. 6 представлены результаты численного расчета спектральной плотности флуктуаций импульса Р при к = 0,1 и Р = 0. Видно, что представленный график носит характер резонансной кривой.

Рис.6. Характерный вид спектральной плотности флуктуации импульса Ор при к = 0,1 и силе Р = 0

СДсо)

80 60 40 20 О

О 0,2 0,4 0,6 0,8 СО

Рис. 7. Характерный вид спектральной плотности флуктуации координаты Оу при к = 0,4 и силе Р = 0

Изображенный на рис. 7 график спектральной плотности флуктуа-ций координаты У при к = 0,4 и Р = 0 помимо резонансной области имеет особенность в области низких частот. Степенная кривая аппроксимации для первых двухсот значений спектральной плотности (ш) описывается выражением

Су (ш) = 0,778ш-0'945, (37)

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

На рис. 8 показана зависимость резонансной частоты от значения параметра к. Видно, что резонансная частота растет с ростом к и может быть аппроксимирована степенной зависимостью

Шрез = 1,157к °'492, (38)

т.е. показатель степени близок к 0,5, как и в случае классического осциллятора. Отметим, что коэффициент пропорциональности в выра-

0,6 0,5 0,4 0,3

0,2 -1-----Г1

0,05 ОД ОД 5 0,2 0,25 к

Рис. 8. Зависимость резонансной частоты осциллятора от параметра к (1) и степенная функция аппроксимации (2)

Gp(га)

J* У / ч ч \ ч

'/ / \ ч \

/ J / \ ч \ Ч

/ / 2 \ V— Ч

/ / // // // ч \ ч

Рис. 9. Резонансные кривые для осциллятора, находящегося под воздействием флуктуирующего коэффициента трения (1), и классического осциллятора (2)

жении (38) близок к 1, что является следствием того, что при расчетах принято М = 1.

На рис. 9 представлены резонансные кривые для классического осциллятора и осциллятора, описываемого выражением (28). В обоих случаях принято, что параметр упругости к = 0,1. Для неклассического осциллятора резонансная кривая получена путем аппроксимации численных результатов полиномом четвертой степени. Видно, что кривая, соответствующая неклассическому случаю, в области средних частот заметно выше резонансной кривой, соответствующей классической теории, а резонансная частота в первом случае несколько больше. Отметим, что суммарная площадь под кривыми в обоих случаях близка к единице.

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

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

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

СПИСОК ЛИТЕРАТУРЫ

1. Бочков Г. Н., Кузовлев Ю. Е. Новое в исследованиях 1/f-шума // Успехи физических наук. - 1983. - Т. 141, вып. 1. - С. 151-176.

2. H o o g e F. N., Kleinpenning T. G. M., V a n d a m m e L. K. J. Experimental studies on 1/f noise // Reports on progress in physics. - 1981. -V.44. - No. 5. - P. 479-532.

3.Пугачев В. С., С и н и ц ы н И. Н. Стохастические дифференциальные системы. - М.: Наука, 1990. - 632 с.

4. Г а р д и н е р К. В. Стохастические методы в естественных науках. - М.: Мир, 1986. - 528 с.

5. Хорстхемке В., Лефевр Р. Индуцированные шумом переходы. -М.: Мир, 1987.-400 с.

6. З у б а р е в Д. Н. Неравновесная статистическая термодинамика. - М.: Наука, 1971.-416 с.

7. Климонтович Ю. Л. Статистическая физика. - М.: Наука, 1982. - 608 с.

8. М о р о з о в А. Н. Необратимые процессы и броуновское движение. - М.: Изд-во МГТУ им. Н.Э. Баумана, 1997. - 332 с.

9. М о р о з о в А. Н. Метод описания немарковских процессов, задаваемых линейным интегральным преобразованием // Вестник МГТУ им. Н.Э. Баумана. Сер. "Естественные науки". - 2004. - № 3. - С. 47-56.

Статья поступила в редакцию 24.12.2007

Андрей Николаевич Морозов родился в 1959 г., окончил МВТУ им. Н.Э. Баумана в 1982 г., д-р физ.-мат. наук, профессор, заведующий кафедрой "Физика" МГТУ им. Н.Э. Баумана. Автор более 100 научных работ в области прецизионных измерений и физической кинетики.

A.N. Morozov (b. 1959) graduated from the Bauman Moscow Higher Technical School in 1982. DSc (Phys.-Math.), professor, head of "Physics" Department of the Bauman Moscow State Technical University. Author of more than 100 publications in the field of high precision measuring systems and physical kinetics theory.

Алексей Владимирович Скрипкин родился в 1983 г., окончил Курский государственный университет в 2005 г. Аспирант кафедры "Физика" МГТУ им. Н.Э. Баумана. Автор 8 работ в области физической кинетики и взаимодействия электромагнитных полей с проводящими средами.

A.V. Skripkin (b. 1983) graduated from Kursk State University in 2005. Postgraduate student of "Physics" Department of the Bauman Moscow State Technical University. Author of 8 publications in the field of physical kinetics theory and interaction of electromagnetic waves with conducting media.