Научная статья на тему 'К вопросу о расчете контактных напряжений и нагрузочной способности эвольвентных зубчатых передач'

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

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

Аннотация научной статьи по физике, автор научной работы — Колосова Е. М., Короткин В. И., Чебаков М. И.

Работа выполнена при частичной финансовой поддержке РФФИ (гранты 06-08-00454, 08-08-00873). Приведены результаты численного моделирования контактного взаимодействия зубьев эвольвентной прямозубой зубчатой передачи при перекосах осей зубчатых колес. Показано отсутствие дополнительного (помимо известного герцевского) влияния приведенного радиуса кривизны контактирующих поверхностей на нагрузочную способность контакта, что опровергает данные, содержащиеся в некоторых литературных источниках. Установлено, что расчет контактной прочности эвольвентных передач по стандартной методике дает заниженные показатели по напряжениям и завышенные по нагрузочной способности в сравнении с результатами, полученными на расчетной модели. Ил. 6. Табл. 4. Библиогр. 6 назв.

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

Текст научной работы на тему «К вопросу о расчете контактных напряжений и нагрузочной способности эвольвентных зубчатых передач»

УДК 621.833: 539.3

К ВОПРОСУ О РАСЧЕТЕ КОНТАКТНЫХ НАПРЯЖЕНИИ И НАГРУЗОЧНОЙ СПОСОБНОСТИ ЭВОЛЬВЕНТНЫХ ЗУБЧАТЫХ ПЕРЕДАЧ

© 2008 г. Е.М. Колосова, В.И. Короткий, М.И. Чебаков

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

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

Предварительные замечания, касающиеся плоской задачи Герца

Как известно, строгое аналитическое решение плоская задача Герца имеет при рассмотрении сжатия упругих цилиндров бесконечной длины с параллельными осями. Тем не менее с минимальной долей условности она используется как основа стандартных инженерных расчетов [3] эвольвентных передач с зубьями конечной длины Ь.

Обозначив q 0 = Fn / Ъ (Н/мм) - погонную нагрузку (называемую интенсивностью нагрузки), распреде-

ленную равномерно по длине Ъ цилиндра (зуба), а н0 -герцевское напряжение (МПа), р - приведенный радиус кривизны контактирующих поверхностей зубьев, запишем для стальных цилиндров (при одинаковом модуле упругости Е = 2 -105 МПа и коэффициенте Пуассона V = 0,3) [4]:

а н о = 186,97 q о/ р

(1)

откуда из условия a H 0 = const следует линейная связь (в данном случае прямая пропорциональность) между q 0 и р , а также (с учетом b =const) между Fn и р .

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

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

Таким образом, с ростом угла а от а 1 до а 2 радиус р возрастет в отношении к р = sin а 2/sin а 1 и, на основании (1), в таком же отношении (при a H0 = = const) можно увеличить нормальную силу Fn , т.е. нагрузочную способность контакта. При этом нагрузочная способность зубчатой передачи, определяемая окружной силой Ft = Fn cos а, увеличится пропорционально sin 2а , т.е. в отношении sin 2а 2 / sin 2а 1.

Перекосы в зацеплении зубчатых передач приводят к неравномерности нагрузок и напряжений вдоль зубьев, при этом длина а площадки контакта может оказаться: а) равной длине зуба (a = b, рис. 1 а, в); б) меньше длины зуба (a < b, рис. 1 г).

б

Рис. 1. Условная схема соотношения длины а площадки контакта и длины Ъ зуба при перекосах осей зубчатых колес

а

в

г

Напряжения при перекосах достигают максимума (а тах ) У одного из торцов и минимума у противоположного, причем для случая а < Ь имеем а тЬ = 0. Схема пятна контакта на рис. 1 а соответствует случаю отсутствия перекоса осей, когда на пятне контакта максимальные напряжения будут равны герцевско-му напряжению а н 0 .

Постановка статической контактной задачи теории упругости

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

Пусть начало прямоугольной системы координат х, у, z находится в точке О первоначального касания, ось х совпадает с биссектрисой угла перекоса осей цилиндров, а перпендикулярная ей ось z находится в плоскости, проходящей через оси цилиндров. Цилиндры имеют соответственно радиусы г и R , а длина их равна Ь . Предполагаем, что торцы цилиндров и их цилиндрические поверхности вне зоны контакта свободны от напряжений, а в зоне контакта отсутствует трение. Цилиндры прижимаются силами Fn, перпендикулярными оси Ох и лежащими в плоскости Ozх, и удерживаются моментами М таким образом, что оси цилиндров сохраняют свое направление в процессе деформации.

Fn

чимся взаимодействием двух четвертей рассматриваемых цилиндров, как это показано на рис. 3.

При конечно-элементном моделировании граничных условий и нагружения в программном пакете ANSYS было задано жесткое закрепление плоской грани В1С1Е1D1, при этом плоская грань В 2 С 2 Е 2 D 2 после приложения нагрузки перемещается поступательно, на гранях A1OD1В1 и А2 OD2В2 заданы условия симметрии, а на линии В 2 С 2 задано отсутствие горизонтального перемещения. Сила Fn приложена к грани В 2 С 2 Е 2 D 2 в направлении, противоположном

направлению оси z . При построении конечно элементной модели взаимодействующих цилиндров в окрестности области контакта использовалось каноническое разбиение на элементы (рис. 4).

С2

En

А

J

r а /

У < O

✓ * R

-Л у-' On

У , ' ч

\

b

Fn

Рис. 2. Схема сжатия силой Fn двух упругих круговых цилиндров

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

Рис. 3. Твердотельная модель взаимодействующих цилиндров

Рис. 4. Конечно-элементная модель

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

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

Результаты расчетов

С целью оценки влияния приведенного радиуса р = ^ / (г + R) кривизны на параметры контакта расчет проводился для двух модификаций прямозубой эвольвентной передачи модуля 3 мм, с числами зубьев

пары 27/34: модификация 1 - а 1 = 20 °, что соответствует р 1=7,72 мм (диаметры моделирующих цилиндров г = 13,85 мм, R = 17,45 мм); модификация 2 -а 2 = 28 °, что соответствует р 2 =10,60 мм (диаметры моделирующих цилиндров г = 19,02 мм, R = 23,95 мм).

(В дальнейшем всем параметрам, соответствующим модификации 1, будем присваивать индекс «1», а модификации 2 - индекс «2»). Перекос осей цилиндров определяли углом у , которому задавали значения 0,001 рад. и 0,0005 рад.

Диапазон усилий Fn1 для расчетов модификации

1 назначали: 500; 1500; 3000; 5000 (Н). Для расчетов модификации 2 диапазон усилий назначали с учетом соотношения Fn2 = крFn1 = 1.373Fn1, т.е. соответственно 686; 2059; 4118; 6860 (Н), где кр есть отношение радиусов кривизны соответственно второй и первой модификаций рассмотренной передачи.

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

1420

1120

821

523

226

77.4 0

1

2 /

L ____—

/

_j

12

Рис. 5. График распределения нормальных (1) и эффективных (2) напряжений по длине цилиндра для случая a < b (расчетный вариант № 3 - табл. 1)

Отметим, что для контроля изложенной здесь схемы расчетов проводились сравнения с результатами расчетов, полученных по формуле плоской задачи Герца (1) в случае отсутствия перекосов осей. Например, при у = 0 , b = 13 мм , р = 10,6 мм и Fn = 2059 Н было получено a Hmax =709 МПа, при

р = 7,72 мм, Fn = 3000 Н - а нтах =1004 МПа, а при р = 7,72 мм, Fn = 5000 Н - а нтах =1322 МПа. Формула Герца (1) дает соответственно а н0=723 МПа, а н0 =1023 МПа и а н0 =1319 МП. Кроме того, расчеты показали, что при увеличении силы Fn в п раз максимальные контактные напряжения увеличиваются, приблизительно, в 4п раз (с погрешностью не выше 3-5 %), формула (1) при этом дает увеличение точно в 4п раз. Изложенное свидетельствует о хорошей точности получаемых результатов.

1670

1400

1140

874

609

345

1

Л 1

2 -

------ —-

2.6

5.2

10.4

х, мм

Рис. 6. График распределения нормальных (1) и эффективных (2) напряжений по длине цилиндра для случая а = Ь (расчетный вариант № 18 - табл. 2)

В табл. 1 представлены результаты расчета по программе «ANSYS» максимальных напряжений для случая а < Ь. Индекс «Н» относится к нормальным напряжениям, индекс «е» - к эффективным. Приведены также результаты расчета интенсивности нагрузки, понимаемой как дА = Fn / а . Для удобства расчетные варианты в таблице пронумерованы.

При сопоставлении соответственных результатов по модификациям 1 и 2 (№ 1-9, 2-10, 3-11, 4-12, 5-13, 6-14, 7-15 и 8-16) видно, что разница в значениях напряжений весьма мала, находясь практически в пределах точности счета. Это свидетельствует о том, что взаимосвязь между Fn и р, как и в задаче Герца (1), близка к прямой пропорциональности.

Из (1) вытекает нелинейная герцевская связь

а н0~ Fnm, в которой показатель степени т = 0,5. Подобная степенная связь может быть приближенно использована и по отношению к а нтах из табл. 1, если положить т « 0,25. Более слабая, чем по Герцу, зависимость напряжения от нагрузки - следствие того, что с изменением последней одновременно меняется длина а площадки контакта (табл. 1), и, следовательно, интенсивность дА нагрузки перестает быть прямо пропорциональной самой нагрузке, а выполняется, приблизительно, пропорциональная зависимость

между величинами дА и F<0'5. Следовательно, интенсивность дА является определяющим фактором, нарушающим герцевские взаимосвязи при перекосах и влияющим на нагрузочную способность контакта.

стн, °е? МПа

стн, Ml 1а

х, мм

Сопоставляя соответственные двум модификациям варианты № 17-21, 18-22, 19-23 и 20-24, приходим к заключению, что и в этом случае, как в рассмотренном выше, сохраняется почти прямая пропорциональность между действующей нагрузкой и приведенным радиусом кривизны.

Кроме того, как и в предыдущем случае, здесь также наблюдается приблизительно пропорциональная зависимость между величинами а н тах и ^0'25. Только теперь она объясняется не изменением длины площадки контакта (которая в сторону увеличения ограничена длиной цилиндра), а перераспределением интенсивности нагрузки и напряжений по длине цилиндра.

Таблица 1

Значения напряжений для случая а < Ь

№ варианта р, мм У, рад •^С^Х Н а, мм qA, Н/мм °ffmax, МПа CTemax, МПа

1 500 5,5 91 885 530

2 0,001 1500 9,4 160 1182 672

3 3000 13,2 227 1420 844

4 7,72 5000 17,0 294 1614 922

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

5 (мод, 1) 500 7,8 64 715 402

6 0,0005 1500 13,7 109 967 520

7 3000 19,2 156 1164 625

8 5000 25,0 200 1343 726

9 686 6,2 111 840 499

10 0,001 2059 11,1 185 1127 673

11 4118 15,8 261 1321 735

12 10,60 6860 21,0 327 1525 839

13 (мод. 2) 686 9,2 75 676 381

14 0,0005 2059 16,4 126 909 519

15 4118 23,4 176 1098 583

16 6860 29,7 231 1297 716

Таблица 2

Значения напряжений для случаяа а = b

№ варианта р, мм У, рад Fn1(Fn2), Н b, мм ÜHmax, МПа CTemax, МПа

17 7,72 (мод. 1) 0,001 3000 13 1430 831

18 5000 1670 1000

19 0,0005 3000 18 1168 655

20 5000 1377 794

21 10,60 (мод. 2) 0,001 4118 13 1354 825

22 6860 1602 1009

23 0,0005 4118 18 1105 644

24 6860 1322 800

Отметим, что расчет по формулам работы [1] дает значительно меньшие напряжения и большие величины длины пятна контакта по сравнению с результатами, представленными в табл. 1. Например, для варианта 12 по формулам работы [1] получаем а н тах = = 1246 МПа, а = 29,2 мм, а для варианта 5 - а н тах = = 637 МПа, а = 11,2 мм.

Обратимся теперь к случаю а = Ь . Расчетные варианты для данного случая легко получаются из вариантов табл. 1, если назначить длину Ь цилиндра равной или меньшей длины а площадки контакта. Выполнив это, к примеру, для вариантов № 3, 4, 7, 8, 11, 12, 15 и 16, представим результаты в табл. 2.

Это перераспределение ведет к тому, что с ростом нагрузки напряжения по длине выравниваются, стремясь в пределе к минимальным (герцевским).

Как показали расчеты, для обоих случаев выявленные и отмеченные выше взаимосвязи между параметрами контакта (кривизны - усилия - напряжения и т.д.) сохранялись для любых значений угла перекоса, вплоть до у = 0,005 рад.

Полезно теперь провести сравнение результатов расчета нагрузочной способности контакта, полученных на модели, и по стандартной методике [3] для обоих рассмотренных случаев. При этом будем полагать отсутствие приработки поверхностей зубьев и принимать во внимание первоначальный (доприрабо-точный) коэффициент Кн р концентрации нагрузки [3].

В соответствии с [3] нормальные контактные напряжения ст Н определяют по формуле

ст Н = ст,

' но^Кнр ,

где ст н 0 - герцевское напряжение (1).

В свою очередь, коэффициент Кнр концентрации

нагрузки при отсутствии скручивания вала, динамических явлений, приработки и т.д. определяется (в принятых в данной статье параметрах) как [3]:

Кн р= 1 + 0,46 2 ус / Fn,

(2)

(3)

Указанное выше сравнение проведем, используя зависимости (2), (3). Тогда при заданных р и Ь запишем с точностью до постоянных коэффициентов:

72

H max

.2

где с - удельная жесткость пары зубьев, которая для модификации 1 нашего примера с1 = 16100 Н/мм2 [3]. Численным моделированием установлено, что с точностью до 3 % можно принять с1 = с 2 = с .

пК нр и ст Н = 1 пКН р , (4)

где Кп и Кнр соответствуют напряжению ст Н , а К'п и К 'н р - напряжению ст н тах , если его рассчитывать по (2).

Обозначая kст = ст нтах / стн , на основании (3), (4) получим:

К = (1п +^ст -<, Кнр = / 1п, К'нр = 1+г / Кп, (5) где t = 0,4Ь2ус.

Сопоставление нагрузочной способности контакта по стандартному расчету и по модельной задаче осуществим с помощью коэффициента

КК = Кп / Кп .

Говоря о варианте а < Ь, отметим, что, согласно существующему регламенту [5], каждой степени точности передачи должна соответствовать минимальная относительная длина площадки контакта, выражаемая как аь = а / Ь (%) и определяемая под легкой нагрузкой на контрольно-обкатных станках. Для весьма грубых степеней точности величина аь должна быть не менее 25 %, для более точных она повышается. При эксплуатационных нагрузках соответствующие реальные цифры будут, естественно, значительно выше, но для полноты картины рассмотрим весь диапазон - от 25 до 100 %.

В табл. 3 представлены выборочные сравнительные результаты для случая а < Ь, а в табл. 4 - для случая а = Ь.

Таблица 3

Сравнительные результаты для случая а < Ь

№ варианта ab, % b, мм t, Н Kw стн, МПа К F n, Н к'нр KF

№ 2 Fn =1500 Н, а=9,4 мм 25 37,6 9105 7,07 1130 1,046 2498 4,64 1,67

50 18,8 2276 2,52 954 1,239 3521 1,65 2,35

70 13,4 1156 1,77 947 1,248 2982 1,39 1,99

100 9,4 569 1,38 998 1,184 2333 1,24 1,56

№ 6 Fn =1500 Н, а=13,7 мм 25 54,8 9670 7,45 961 1,006 1640 6,90 1,09

50 27,4 2417 2,61 804 1,203 3249 1,74 2,17

70 19,6 1237 1,82 794 1,218 2823 1,44 1,88

100 13,7 604 1,40 833 1,161 2231 1,27 1,49

№ 10 Fn =2059 Н, а=11,1мм 25 44,4 12696 7,17 1047 1,076 4400 3,89 2,14

50 22,2 3174 2,54 881 1,279 5389 1,59 2,62

70 15,9 1628 1,79 874 1,289 4503 1,36 2,19

100 11,1 793 1,39 922 1,222 3468 1,23 1,68

№ 14 Fn =2059 Н, а=16,4 мм 25 65,6 13857 7,73 894 1,017 2598 6,33 1,26

50 32,8 3464 2,68 745 1,220 4758 1,73 2,31

70 23,4 1763 1,86 734 1,238 4099 1,43 1,99

100 16,4 866 1,42 766 1,187 3253 1,27 1,58

Таблица 4

Сравнительные результаты для случая a = Ь

№ варианта Ь,мм t, Н кщ uH, МПа К Fт Н к'нр KF

17 13 1088 1,36 1192 1,200 4799 1,23 1,60

18 1,22 1457 1,152 6995 1,16 1,40

19 18 1043 1,35 1009 1,158 4375 1,24 1,46

20 1,21 1233 1,117 6491 1,16 1,30

21 13 1088 1,26 1147 1,181 6168 1,18 1,50

22 1,16 1420 1,128 9024 1,12 1,32

23 18 1043 1,25 971 1,138 5644 1,18 1,37

24 1,15 1202 1,100 8520 1,12 1,24

По результатам из табл. 3 и 4 отчетливо видно, что нагрузочная способность контакта, рассчитываемая по стандартной методике [3], существенно выше, чем рассчитанная по модели. Есть основания полагать, что эта разница возрастет, если сравнительную оценку производить не по нормальным, а по эффективным напряжениям, которые для пространственных контактных задач считаются в инженерной практике критериальными [6]. Это следует из того, что при использовании стандартной методики, основанной на плоской задаче Герца, отношение эффективного напряжения к нормальному при V = 0,3 равно 0,4, в то время как по данным табл. 1 и 2 отношение ст етах/ стнтах колеблется в диапазоне 0,53...0,63. Данное превышение логично, если учесть, что вблизи торца резко уменьшается (или вовсе исчезает) одно из главных напряжений.

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

Что же касается используемого в стандартном расчете коэффициента Кнр (КНр) концентрации

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

Выводы

1. При сжатии двух упругих круговых цилиндров с пересекающимися осями (модельная задача) связь между сжимающим усилием и приведенным радиусом

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

2. Определяющим параметром, нарушающим гер-цевскую связь усилия с напряжением при перекосах осей цилиндров и существенно влияющим на нагрузочную способность контакта, является изменяющаяся под нагрузкой так называемая ее интенсивность дА = Fn / а.

3. Расчёт контактной прочности эвольвентных зубчатых передач по ГОСТу дает более низкие контактные напряжения и более высокую нагрузочную способность контакта, чем расчет по модельной задаче.

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

зависящий.

Работа выполнена при частичной финансовой поддержке РФФИ (гранты 06-08-00454, 08-08-00873).

Литература

1. Гришин С.А. Контактное взаимодействие упругих цилиндров при перекосах осей // Респ. сб. «Теоретическая и прикладная механика». Вып. 19. Харьков,. 1988. С. 32-39.

2. Журавлев Г.А. Оценка применимости решения Герца в задачах о контакте зубчатых колес // Техника машиностроения. 2001. № 2. С. 82-90.

3. ГОСТ 21354-87. Передачи зубчатые цилиндрические звольвентные. Расчет на прочность. М., 1988.

4. Кудрявцев В.Н., Державец Ю.А., Глухарев Е.Г. Конструкции и расчет зубчатых редукторов. Справочное пособие. Л., 1971.

5. Г0СТ 1643-81. Передачи зубчатые цилиндрические. Допуски. М., 1981.

6. Ковальский Б.С. Расчет деталей на местное сжатие. Харьков, 1967.

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

НИИ механики и прикладной математики

им. И.И. Воровича Южного федерального университета 21 апреля 2008 г.

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