Научная статья на тему 'Нестационарные колебания растущей круговой цилиндрической оболочки'

Нестационарные колебания растущей круговой цилиндрической оболочки Текст научной статьи по специальности «Физика»

CC BY
122
44
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
РАСТУЩАЯ ЦИЛИНДРИЧЕСКАЯ ОБОЛОЧКА / ВЫНУЖДЕННЫЕ КОЛЕБАНИЯ / АНАЛИТИЧЕСКИЕ РЕШЕНИЯ / СПЕКТРАЛЬНЫЕ РАЗЛОЖЕНИЯ / СОБСТВЕННЫЕ ФУНКЦИИ / GROWING CYLINDRICAL SHELL / NON-STATIONARY VIBRATIONS / ANALYTICAL SOLUTIONS / SPECTRAL DECOMPOSITION / EIGENFUNCTIONS

Аннотация научной статьи по физике, автор научной работы — Барышев А. А., Лычев С. А., Манжиров А. В.

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

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Non-Stationary Vibration of Growth Circular Cylindrical Shell

Small forced vibrations of growing cylindrical shell fixed on circular boundaries is studied in the framework of Kirchhoff-Love shell theory. The process of the accretion are characterized by the continuous adherence of material particles to its facial surface. Since the shell bends during the accretion, its stressed-strained state depends not only on loading, but also on the history of the process of accretion, i.e. the schedule of accretion. A schedule of elementary type, when during every infinitesimal time interval the particles of adhered material constitute the layer of constant infinitesimal thickness, is considered.

Текст научной работы на тему «Нестационарные колебания растущей круговой цилиндрической оболочки»

МЕХАНИКА

УДК 539.3

НЕСТАЦИОНАРНЫЕ КОЛЕБАНИЯ РАСТУЩЕЙ КРУГОВОЙ ЦИЛИНДРИЧЕСКОЙ ОБОЛОЧКИ

А. А. Барышев, С. А. Лычев , А. В. Манжиров

Саратовский государственный университет E-mail: BaryshevAA@gmail.com

* Институт проблем механики им. А. Ю. Ишлинского РАН, Москва E-mail: lychevsa@mail.ru, manzh@inbox.ru

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

Ключевые слова: растущая цилиндрическая оболочка, вынужденные колебания, аналитические решения, спектральные разложения, собственные функции.

Non-Stationary Vibration of Growth Circular Cylindrical Shell A. A. Baryshev, S. A. Lychev, A. V. Manzhirov

Small forced vibrations of growing cylindrical shell fixed on circular boundaries is studied in the framework of Kirchhoff-Love shell theory. The process of the accretion are characterized by the continuous adherence of material particles to its facial surface. Since the shell bends during the accretion, its stressed-strained state depends not only on loading, but also on the history of the process of accretion, i.e. the schedule of accretion. A schedule of elementary type, when during every infinitesimal time interval the particles of adhered material constitute the layer of constant infinitesimal thickness, is considered.

Key words: growing cylindrical shell, non-stationary vibrations, analytical solutions, spectral decomposition, eigenfunctions.

ВВЕДЕНИЕ

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

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

© Барышев А. А., Лычев С. А., Манжиров А. В, 2012

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

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

1. ПОСТАНОВКА ЗАДАЧИ

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

1 + V д2ь

/ д2 1 - ид2

НдО? + ~

и + Б-

в

дв2] ” ' _ 2 1- ид2

V дт д2 и

дадв + Я да Р дЬ2 ‘

+ Б

д2

да2 + дв2

1 + и д2и 2 дадв и ди Б дь {Б 2 2

БЯда + Ядв + ( я2+БУ у

Б дт д2ь

Ь + Ядв = РН~Щ2

и,+рнд22 = <э.

Здесь а и в — координаты на срединной поверхности, отсчитываемые вдоль образующей и в окружном направлении (рис. 1), р — плотность массы, ЕН ЕН3

Б =

1 — и2 ’

Б =

12(1 - и2)

— тангенциальная

и изгибная жесткости, Е и и — модуль Юнга и коэффициент Пуассона материала, и, ь и т — тангенциальные смещения и прогиб точек срединной поверхности, Ь — время, = <^(Ь) — интенсивность распределенной по лицевой поверхности поперечной нагрузки.

Будем считать, что оболочка замкнута, а ее края жестко заделаны. В этом случае краевые условия на основаниях формулируются в виде

Рис. 1. Расчетная схема оболочки: а — оболочка постоянного состава, б — растущая оболочка

дт

да

а=0,2Ь

= 0,

и\а=0,2Ь = ь\ а=0,2Ь = т\а=0,2Ь =

к которым следует добавить условия периодичности:

и(в + 2пЯ) = и (в), ь(в + 2пЯ) = ь(в), т(в + 2пЯ) = т(в).

(1)

(2)

Начальные условия определяют в начальный момент времени тангенциальные смещения и0,ь0 и прогиб т0 точек срединной поверхности, а также их начальные скорости и0,г}0,т0, а именно:

ди дь дт

,0 = 0 \ ,0 = 0 ,0 т0 = 0 т\ ~дЬ = Щ, t=0 д = $0, t=0 ~дЬ

t=0

= 1^0 ■

(3)

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

/ д2 1 - ид2

уда2 + 2 дв2

_ 1 + и д2ь ^ и дт ,...

и+Б~ шв + Бяэа =рНи-

(4)

1 + V д2ь . - - „ „ . . -

В----1-- „ „ „ + В ( ---— 7^7 + ^7 ) V + — — = рН) ,

2 дадв \ 2 да2

вV дй В ді) / В

яда + Ядр + V Я2

1 - V д2 д2

дв2

+ БУ2V2 ) Ъ + рНЪ = 2.

В дЪ

яде

(5)

(6)

В (4)-(6) точкой обозначена производная по времени. Отметим, что уравнения растущих оболочек имеют третий порядок производных по времени, а параметр Н, определяющий толщину и входящий в коэффициенты, меняется со временем Н = Н(£).

Считается, что скорость присоединения слоя совпадает со скоростью движения лицевых поверхностей оболочки. В противном случае, в уравнения (4)-(6) необходимо добавить слагаемые, характеризующие тормозящую силу, возникающую вследствие учета силы инерции присоединяемого материала в процессе роста.

Краевые условия сохраняют вид (1)-(2), а начальные данные в отличие от (3) содержат начальные ускорения йо, ^о, Шо точек срединной поверхости оболочки и записываются в следующей форме:

и|*=о = йо, у|*=о = Уо, -ш1г=о = Шо,

дй д) дъ д 2й 2 д 2ъ

д = йу0, і=0 д = У0, і=0 ~ді = Ъ0, і=0 ді2 = й0, і=0 ді2 = У0, і=0 ді2

і=0

= ъо.

Для построения решения сформулированных задач целесообразно перейти к безразмерным величинам, а именно безразмерным координатам — х = а/Ь — 1, р = @/Я; безразмерному времени —

1 I Е „ „

т = -----— £; безразмерным тангенциальным смещениям и прогибу — й = й/Я, У = у/Я,

НЯ

Я\ р(1 — V2)

ъ = ъ/Я. Будем использовать безразмерные параметры 5 = Я/Ь, у = __

V12Ь2

~ д2 д2 тор Лапласа V2 = 7—77 + 5—2—^.

дх2 д^2

Введем оператор £, характеризующий упругую реакцию оболочки на перемещения:

и безразмерный опера-

(— (52

=

д2 , 1-у д2

+ 2 др2

дх2

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

5

д2

V

2 дхдр

V5 дХ

1+У 5 д

2 дхдр

— ( 1-V 52 + _д2_

1 2 5 дх2 + др2

д др

—V5 д ^

дх

д

др

1 + Н2 'V2 'V2)

область определения которого задается оператором краевых условий:

В [и] = ( й|х=-1 ,г;|х=-1,Ъ;|х=-1,

х = 1

(7)

В формуле (7) используется вектор смещения с безразмерными компонентами: И = (й, V, Ъ). Определим также вектор начального смещения И0 = (й0,)0,ъ0) /Я, вектор начальной скорости

Уо =

р(1 — V2) Е

(й0,У0,Ъ0), вектор начального ускорения W0 =

р(1 — V2) Е

Я (й0,у0, Ъ0), вектор

1 — V

безразмерных интенсивностей нагрузки Q = ——— Я (0,0,2) и вектор скоростей интенсивности

ЕН

1 — V2 / -\

у = ~ттЯ (°,0’2) •

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

£ [И]+ И = Q, В [И]=0, И|т=0 = И0, И|т=0 = У0, (8)

динамическая задача о вынужденных колебаниях растущей оболочки (задача 2):

£ |И1 + И = , В[И]=0, И|т=0 = И0, И|т=0 = У0, И|т=0 = Wо. (9)

Соотношения (8), (9) определяют математические формулировки рассматриваемых в работе задач о вынужденных колебаниях нерастущей и растущей цилиндрической оболочки.

2. ПОСТРОЕНИЕ РЕШЕНИЯ В ФОРМЕ СПЕКТРАЛЬНЫХ РАЗЛОЖЕНИЙ

Решение начально-краевых задач (8), (9) будем искать в классе интегрируемых с квадратом вектор-функций, определенных в области [—1; 1] х [0; 2п] в форме разложения по собственным функциям оператора £ [3]. Эти функции являются нетривиальными решениями задачи Штурма-Лиувилля:

L [U] = AU, B [U] = 0.

(10)

В силу самосопряженности оператора L его собственные значения являются вещественными, а собственные функции — ортогональны и образуют полную в пространстве L2 систему функций.

Условиям периодичности по окружной координате p для любого целого m удовлетворяют собственные функции, образующие два семейства: четные по нормальным перемещениям — индекс «е» (even) и нечетные — индекс «о» (odd).

Ue =

v*/ em

(

V

Unm COS mp

v°nm sin mp

\

w“ COS

(

U° = em

\

V

u°m sin mp v°°m cos mp

mp J

w" sin

В дальнейшем индексы е, о будем опускать, поскольку процедуры нахождения собственных функций идентичны. Формулы, определяющие собственные функции ипт = (йпт,Упт,Шпт), в явном виде приведены в приложении.

го

Решение для оболочки постоянного состава будем искать в виде И(х,р,т) = ^ 'фпт(т) х

п,т=о

го

х Ипт(х, р), а для растущей оболочки — в форме И(х,р,т) = ^ фпт(т)Ипт(х, р,т).

п,т=о

Заметим, что в случае вынужденных колебаний растущей оболочки решение не предполагает отделения пространственных переменных и времени.

Координатные функции ,фпт(т) задачи 1 находятся из бесконечной последовательности несвязанных задач Коши:

Фы + AkiФы = qki, Фы(0) = (Uq, Uki), фкь(0) = (Vq, Ukl}.

(11)

В задаче о растущей оболочке для определения координатных функций фпт (т) приходим к бесконечной последовательности связанных задач Коши:

+ Akl Фы + УУ

°,m=Q

dUnm тт \ „ ,, /d2Unm тт

QT , Ukl) + 3Фет\ дТ2 , Uk^ +

/ d3U dU

, ^ / °m | \ WKJnm TT

+ Ynm { QT3 + A°m , U k

= Qkh

фЫ (0) = (Uq, U

kl}

т=Q

n,m=0

д U°

дт

, Uk

т=Q

= (Vq, Ukl}

т=0

ф’1, (0)+ £ (2ф°„,(0) dUU°m + ф°т (0)d2U°m, Ukl

n,m=0

т=0

= (Wq, Ukl}

т=0

(12)

(13)

(14)

гч | т пт \ ) с\ о

дт дт2

В (11)-(14) угловыми скобками обозначено скалярное произведение функций в пространстве £2, т. е.

1 2п

(U°m, U°m} — II (unm + vnm + Wnm)

-1 0

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

а яы = (Q,ик1 } ш = ^, и^}.

Решение записанных задач может быть получено методом усечения, а численное интергрирование проведено, например, методом Рунге-Кутта.

го

3. АНАЛИЗ ЧИСЛЕННОГО РЕШЕНИЯ

В качестве примера проведен анализ напряженно-деформированного состояния круговой цилиндрической оболочки, выполненной из алюминия. Геометрические размеры и механические свойства материала оболочки приняты следующими: Я = 0.01 м, начальная толщина Л,0 = 10-4 м, Ь = 0.01 м; Е = 68 ГПа, V=0.36, р=2700 кг/м3, а также безразмерная интенсивность внешней нагрузки до = 0.01.

Будем считать, что присоединение материала происходит с постоянной скоростью и толщина оболочки изменяется по закону И,(Ь) = Л,0 + ьдгЬ.

В работе рассмотрен важный с практической точки зрения случай, когда к началу наращивания

1- V2

оболочка изогнута равномерно распределенной нагрузкой постоянной интенсивности д0 =

ЕН

-ЯЯо

В момент начала наращивания нагрузка снимается, и оболочка приходит в движение. Тогда векторы

и0, У0, Wo примут следующий вид: и0 = Е

Яп

п,т=0

-ие

^ пт

(х,р), Уо = 0, Wo = Е

п,т=0

ЯптиПт (х,^

1 2п

где дпт = оо / / Ыпш (х, ф) йрйх.

-1 о

На рис. 2 изображены поверхности, иллюстрирующие собственные функции оболочки при т = 1 и различных значениях п. Так, на рис. 2, б изображена функция, отражающая преимущественно нормальную форму, г — продольную форму, а д — крутильную.

а б в

Рис. 2. Собственные формы колебаний цилиндрической оболочки при п

п = 44 (г), п = 45 (д)

д

0 (а), п = 6 (б), п =18 (в),

Численное определение собственных значений и собственных функций краевой задачи для оператора С имеет значительные технические трудности, связанные с возможным наличием кратных собственных значений и кратных корней алгебраического уравнения (16). В связи с этим контроль получаемого разложения возможно осуществить путем разложения в ряд по собственным функциям ипт вектор-функции (0,0,1). Такие разложения представлены на рис. 3, а-в при т =1 и с различным количеством слагаемых, причем последнее слагаемое имеет номер п0.

\.2[

0.8

0.6

0.4

0.2

1.0 -1.0 -0.5

0.5

б

Рис. 3. Тестовое разложение вектор-функции (0,0,1) при по = 15 (а), по = 30 (б), по = 60 (в)

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

ос

г

а

в

2000

1000

-1000

-2000

0 С 0і

Ю)(

10

ф?20

1000

500

-500

-1000

0:

х'ср

.115

)}Щ20

б

Рис. 4. Координатные функции оболочки при т =1 и п = 3 (а), п = 5 (б)

На рис. 5 изображены изменения во времени нормального смещения точки срединной поверхности оболочки с безразмерными координатами (0.6,0). Они показывают, что присоединение материала к оболочке уменьшает амплитуду колебаний.

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

Рис. 5. Осцилограммы колебаний растущей оболочки и оболочки постоянного состава

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

Работа выполнена при финансовой поддержке РФФИ (проекты 11 -01-00669-а, 12-08-01119-а, 12-08-01260-а).

Библиографический список

1. Биргер И. А., Пановко Я. Г. Прочность. Устойчивость. Колебания : справочник : в 3 т. М. : Машиностроение, 1968. Т. 3. 567 с.

2. Лычев С. А., Лычева Т. Н., Манжиров А. В. Неста-ПРИЛОЖЕНИЕ

ционарные колебания растущей круглой пластины // Изв. РАН. МТТ. 2011. № 2. С. 199-208.

3. Коллац Л. Задачи на собственные значения (с техническими приложениями). М. : Наука, 1968. 504 с.

Приведем явные формулы собственных функций ипт краевой задачи (10) и последовательность их определения. Индексы п, т в записываемых далее формулах будем опускать. Определим мат-

/Лів^іХ А2 в^2 х ......... Л8в^х\

рицу фундаментальных решений М =

В1 е^1 х \Сі е^1 х

В2 е^2х С2е^2х

В8 е^8 С8е^х/

и = МБ.

и вектор постоянных

(15)

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

Лк = трк6 (р - (і + Н2 (рк - 6-2т2)2 - А^

а

Bk = ( 52мк - -1—^m2 + Л ) (1 + h2(мк - 5 2m2)2 - Л) - (v5Mk)2,

Ck = m

1 + V с2 2 Лт2 2 1 - 2

—2— Мк - ( 5 мк

m2 + Л

причем Мк являются корнями характеристического уравнения:

boM8 + bi м6 + b2M4 + Ьз м2 + &4 = О,

(16)

коэффициенты которого определяются равенствами 60 = Л2(1-V)/2,6і = 52Л2(А(3—V)—4т2(1-V))/2, 62 = (Л2(2А2 — 3т2А(3 — V) + 6т4(1 — V)) + 54(1 — V)(1 — V2 — А))/2, 63 = 5-2(-2т2А(2Л2А —

— 54(1 — V)) — 3т4Л2А(—3 + V) — 4т6Л2(1 — V) + 54А(3 — А(3 — V) — V — 2v2))/2, 64 = 5-4(—т6Л2 + + т2 54 А + т4Л2А — 54 (—1 + А)А)(2А — т2(1 — v))/2.

Собственные значения А определяются численно из трансцендентного уравнения, получающегося при действии оператора В на решения (15). Его можно записать следующим выражением:

3. - e Ai 2 L - e 2 A2 . . . A8e-L8

Bi e ■R 2 L - e 2 B2 ... B8e-L8

3. - e Ci 2 L - e C2 ... C8e-L8

CiMie-Ll C2M2e-L2 ... . . . C8M8e-L8

AieLR A2eL2 ... ... A8 eL8

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Bi eLl B2 eL2 ... . . . B8eL8

Ci eLl C2 eL2 ... .. . C8eL8

CiMieLl C2M2eL2 . .. ... C8M8eL8

= О.

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

2втЬ (^ + ^)'

|U||2 = (U, U) = п ( DiDj (AiAj + BiBj + CiCj)

i,j = i

Mi + Mj

УДК 519.68:[5/6+3]; 004.94; 544.015.4; 544.022.822

ПЕРКОЛЯЦИЯ СФЕР В КОНТИНУУМЕ

М. М. Бузмакова

Астраханский государственный университет E-mail: mariya_nazarova@mail.ru

Предложена модель континуальной перколяции жестких сфер с проницаемыми оболочками, которая описывает фазовый переход золь-гель. Сферы имеют жесткие части радиусом г, которые не могут перекрываться друг с другом, и проницаемые оболочки шириной d, которые могут перекрываться. Такие сферы одинакового размера случайным образом помещаются в куб с линейным размером L. Вероятность возникновения связи меж-дусферами пропорциональна объему перекрытия проницаемых оболочек. Если связь между сферами возникает, то сферы принадлежат одному кластеру. В задаче ищется перколяционный кластер, т. е. кластер, соединяющий нижний и верхний грани куба. Доля заполнения куба сферами, при которой вероятность возникновения перколяционного кластера равна 0.5, называется порогом перколяции. Порог перколяции соответствует точке геля. Получена зависимость значения порога перколяции оттол-щины проницаемой оболочки.

Ключевые слова: компьютерное моделирование, перколяция, фазовый переход золь-гель.

Percolation of Spheres in Continuum M. M. Buzmakova

The model of the continuum percolation of hard spheres with permeable shells, which describes phase transition sol-gel,has been investigate. Spheres have hard parts in radii r, which can’t be blocked with each other, and permeable shells in width d, which can be blocked. Such spheres of the equal size have been randomly packing in the cub with linear size L. The probability of joining the spheres in a cluster is proportional to the volume of overlapping of permeable shells. Spheres belong to a cluster, if a communication between spheres arises. The percolation cluster is the cluster connecting bottom and top sides of the cube. The packing fraction, at which probability of occurrence of the percolation cluster is 0.5, is called as the percolation threshold. The percolation threshold corresponds to the gel point. The dependency of the percolation threshold of the hard spheres with permeable shells from a thickness of the shell has been obtained.

Key words: computer modeling, percolation, sol-gel phase transition.

© Бузмакова М. М., 2012

i Надоели баннеры? Вы всегда можете отключить рекламу.