Полный текст

Введение Конструктивные элементы из эллиптических цилиндров являются неотъемлемой частью многих техносферных объектов промышленного и гражданского назначения. Определение прочностных параметров эллиптических цилиндров, особенно в геометрически нелинейной постановке, в настоящее время осуществляется посредством численных методов расчета [1-7], наиболее распространенным из которых является метод конечных элементов (МКЭ) [8-14]. При разработке конечно-элементных алгоритмов расчета эллиптических цилиндров исследователям неизбежно приходится решать вопрос о реализации интерполяционной процедуры, являющейся по своей сути фундаментальным аспектом МКЭ. Если в качестве глобальных координат при расчете эллиптических цилиндров выбрать «естественные» криволинейные координаты, связанные с параметрами срединной поверхности оболочки данного типа, то необходимо использовать векторную форму интерполяционной процедуры [15; 16]. В изложенной работе представлен алгоритм расчета прочности эллиптического цилиндра в геометрически нелинейной постановке с учетом деформаций сдвига при реализации векторной формы интерполяции искомых величин. 1. Материалы и методы Положение точки срединной поверхности эллиптического цилиндра можно задать радиус-вектором (1) где - орты декартовой системы координат; - полярный угол эллипса, являющегося поперечным сечением эллиптического цилиндра плоскостью, перпендикулярной продольной оси ОХ цилиндра. Полярный радиус поперечного сечения эллиптического цилиндра определяется следующей зависимостью Выражения для ортов локального базиса точки срединной поверхности эллиптического цилиндра могут быть получены с использованием формулы дифференциальной геометрии (2) где - длина дуги эллипса поперечного сечения эллиптического цилиндра с параметрами и b. Из (2) можно получить прямую и обратную матричные зависимости (3) где Положение точки расположенной в произвольном слое эллиптического цилиндра, находящегося на расстоянии от срединной поверхности, в недеформированном состоянии определяется радиус-вектором (4) После серии последовательных шагов нагружения точки и переместятся в положения и , которые могут быть определены соответствующими радиус-векторами (5) где и - суммарные векторы перемещения точек и за j шагов нагружения. После завершения (j+1)-го шага нагружения точки M и перейдут в точки M * и M *t , положения которых могут быть определены посредством радиус-векторов (6) где и - векторы шаговых перемещений точек M и . Входящие в (5) и (6) суммарный и шаговый векторы перемещения точки определяются по формулам (7) где и - суммарный и шаговый векторы углов наклона нормали (8) где индексы a последовательно принимают значения 1 и 2. Выполняя операцию дифференцирования (5) и (6), можно получить выражения для векторов базиса в точке произвольного слоя эллиптического цилиндра после j шагов нагружения и на (j+1)-м шаге (9) где индекс принимает значения x и j. Приращение деформаций в произвольном слое эллиптического цилиндра на (j+1)-м шаге нагружения можно получить, используя соотношение механики сплошных сред [17] (10) где Поверхность эллиптического цилиндра моделируется совокупностью четырехузловых элементов дискретизации, ориентированных вдоль линий главных кривизн. Узловыми варьируемыми параметрами дискретного элемента являются компоненты вектора шагового перемещения точки срединной поверхности, первые производные этих компонент и компоненты шагового вектора углов поворота нормали. Столбцы искомых узловых неизвестных в глобальной и локальной системах координат могут быть представлены в следующем матричном виде (11) (12) где Здесь под понимается компонента вектора шагового перемещения , и . (13) В настоящей статье при построении интерполяционных выражений реализованы два варианта. В первом варианте применялся общепринятый в МКЭ подход [18-25], при котором каждая компонента вектора шагового перемещения и каждая компонента шагового вектора углов поворота нормали интерполируются через узловые значения своей же компоненты (14) Входящие в (14) функции формы и представляют собой произведения полиномов Эрмита третьей степени и билинейные функции локальных координат соответственно. Использование стандартной интерполяционной процедуры (14) оправдано при решении задач в прямоугольной декартовой системе координат. Однако при исследовании НДС тел, имеющих криволинейные очертания (к которым относятся оболочки вообще и эллиптические цилиндры в частности) целесообразнее применять криволинейные системы координат (цилиндрические, сферические, тороидальные), специально разработанные для анализа объектов с криволинейными границами. На криволинейной поверхности локальный базис постоянно изменяется при переходе от точки к точке, поэтому интерполяционные выражения должны содержать параметры используемой системы криволинейных координат. Такому критерию удовлетворяет векторная форма интерполяционной процедуры, которая и была реализована во втором варианте расчёта. Узловыми варьируемыми параметрами при этом выбираются непосредственно векторы шагового перемещения, первые производные этих векторов и шаговые векторы углов поворота нормали (15) (16) (17) Применяя интерполяционные выражения (14) к векторным узловым неизвестным, можно записать следующие интерполяционные зависимости (18) (19) где - матрица перехода от (16) к (15). Входящие в (18) и (19) узловые векторы шагового перемещения, а также первые производные этих векторов по криволинейным координатам и , шаговые векторы углов поворота нормали могут быть представлены компонентами узловых векторов локального базиса (20) (21) где верхний индекс обозначает узлы конечного элемента Соотношения (18) и (19) с учётом (20), (21) запишутся следующим образом (22) (23) где Элементами квазидиагональных матриц и являются соответствующие матрицы-строки, содержащие узловые векторы локального базиса (24) Матрицы (24) на основании (2) могут быть выражены через орты декартовой системы координат (25) С учётом (3) матричные выражения (25) примут вид (26) Представляя вектор шагового перемещения w и шаговый вектор углов поворота нормали компонентами локального базиса точки внутренней области конечного элемента и учитывая (26), интерполяционные выражения (22), (23) могут быть записаны в виде (27) (28) где Из сопоставительного анализа левых и правых частей (27) и (28) можно получить интерполяционные выражения для компонент вектора шагового перемещения и компонент шагового вектора углов поворота нормали (29) Интерполяционные зависимости в варианте стандартной для МКЭ (14) процедуры или в варианте новой (29) векторной формы интерполяции искомых величин позволяют скомпоновать следующее матричное соотношение (30) где Матрица жесткости и столбец узловых усилий конечного элемента на (j+1)-м шаге нагружения компонуются путём минимизации функционала Лагранжа согласно стандартной для МКЭ в форме метода перемещений процедуре [26-35]. 2. Пример расчета 1 Была решена задача по исследованию НДС цилиндрической оболочки, загруженной в середине пролёта сосредоточенной силой P и имеющей на диаметрально противоположной стороне шарнирную опору, препятствующую вертикальному смещению. По торцам оболочка имела скользящую заделку. Вследствие наличия плоскостей симметрии была рассчитана 1/4 часть цилиндра. Расчетная схема конструкции представлена на рис. 1. Были приняты следующие исходные данные: длина образующей Первоначально параметры эллипса поперечного сечения цилиндрической оболочки принимались одинаковыми т.е. цилиндр был круговым. Расчёты были выполнены в двух вариантах: в первом варианте при формировании матрицы жесткости и столбца узловых усилий конечного элемента на шаге нагружения была использована стандартная для МКЭ интерполяционная процедура (14); во втором варианте для интерполяции компонент вектора шагового перемещения и шагового вектора углов поворота нормали была применена разработанная векторная форма интерполяции искомых величин (29). Результаты повариантного расчёта кругового цилиндра представлены в табл. 1, в которой приведены численные значения нормальных и сдвиговых напряжений на внутренней и наружной поверхностях оболочки на торце цилиндра в точках и, а также значение прогиба в точке приложения силы в зависимости от густоты сетки дискретизации оболочки при фиксированном числе шагов нагружения, равном 30. Рис. 1. Расчетная схема эллиптического цилиндра, загруженного сосредоточенной силой P Fig. 1. Design diagram of an elliptical cylinder loaded with a concentrated force P Анализ данных табл. 1 показывает, что в обоих вариантах расчета кругового цилиндра наблюдается устойчивая сходимость вычислительного процесса по мере сгущения сетки дискретизации. Контролируемые параметры НДС кругового цилиндра в точках А и В (которые должны быть одинаковыми, исходя из симметрии расчетной схемы) практически совпадают и в первом, и во втором вариантах расчета. Различия между вариантами расчета пренебрежимо малы. В табл. 2 представлены повариантные значения напряжений и прогиба под сосредоточенной силой в концевом и срединном сечениях кругового цилиндра в зависимости от количества шагов последовательного нагружения при фиксированной сетке дискретизации рассчитываемого фрагмента оболочки 25´25 узлов. Анализ данных табл. 2 позволил сделать вывод об устойчивой сходимости вычислительного процесса с увеличением количества шагов нагружения в обоих вариантах расчета. Также можно отметить, что при расчете кругового цилиндра в обоих вариантах расчета наблюдаются практическое совпадение значений напряжений в точках А, В и Р, С соответственно, чего и следовало ожидать, учитывая наличие вертикальной плоскости симметрии расчетной схемы. Результаты расчётов позволяют сделать вывод, что при анализе НДС кругового цилиндра в геометрически нелинейной постановке оба варианта интерполяционной процедуры дают практически одинаковые значения параметров НДС оболочки. 3. Пример расчета 2 Рассмотрим задачу исследования НДС эллиптического цилиндра, расчётная схема которого совпадает с круговым цилиндром из примера 1. Соотношение параметров эллипса поперечного сечения эллиптического цилиндра примем равным Результаты повариантных расчётов представлены в табл. 3, структура которой идентична структуре табл. 1. Анализируя данные табл. 3, можно отметить, что и в первом, и во втором вариантах расчета наблюдается сходимость вычислительного процесса со сгущением сетки дискретизации. Однако между вариантами расчета имеются принципиальные и весьма существенные различия. Так, в первом варианте расчета эллиптического цилиндра значения напряжений в точке В примерно в два раза меньше, чем в точке А, в то время как при расчете кругового цилиндра значения напряжений в данных точках практически совпадали. Причем сгущение сетки дискретизации не уменьшает погрешности вычислений по первому варианту. Значения прогиба под сосредоточенной силой в первом варианте также оказались заниженными в 1,67..1,61 раза по сравнению со вторым вариантом расчета. Во втором варианте расчета можно отметить полное совпадение значений напряжений, вычисленных в точках А и В, что и должно быть при выбранной расчетной схеме эллиптического цилиндра. В табл. 4, структура которой совпадает с таковой табл. 2, можно отметить сходимость повариантных расчетов с увеличением количества шагов нагружения. Но, как и в табл. 3, можно констатировать существенные погрешности вычислений напряжений и прогиба в первом варианте расчета, как в концевом, так и в срединном сечениях эллиптического цилиндра. Таблица 1 / Table 1 Численные значения нормальных и сдвиговых напряжений на внутренней и наружней поверхностях кругового цилиндра Numerical values of normal and shear stresses on the inner in and outer n surfaces of a circular cylinder Координаты точек x, j (м, рад) Напряжения, МПа, прогиб, см Вариант интерполяционной процедуры Скалярная Векторная Сетка узлов 21´21 24´24 26´26 28´28 21´21 24´24 26´26 28´28 13,30 12,96 12,78 12,61 13,16 12,85 12,68 12,53 53,63 53,35 55,22 53,10 53,23 53,04 52,95 52,87 -6,85 -7,64 -8,17 -8,70 -6,80 -7,60 -8,13 -8,67 -12,87 -12,57 -12,40 -12,25 -12,74 -12,47 -12,32 -12,18 -56,88 -56,97 -57,07 -57,17 -56,44 -56,48 -56,83 -56,96 -6,85 -7,64 -8,17 -8,70 -6,80 -7,60 -8,13 -8,66 13,29 12,96 12,77 12,61 13,16 12,85 12,68 12,53 53,61 53,34 53,21 53,10 53,23 53,04 52,95 52,87 -6,85 -7,64 -8,17 -8,70 -6,80 -7,60 -8,13 -8,67 -12,88 -12,57 -12,41 -12,25 -12,74 -12,47 -12,32 -12,18 -56,84 -56,98 -57,08 -57,17 -56,44 -56,68 -56,83 -56,96 -6,85 -7,64 -8,17 -8,70 -6,80 -7,60 -8,13 -8,66 -0,5717 -0,5748 -0,5763 -0,5776 -0,5738 -0,5764 -0,5777 -0,5787 Таблица 2 / Table2 Численные значения нормальных и сдвиговых напряжений на внутренней и наружней поверхностях кругового цилиндра Numerical values of normal and shear stresses on the inner in and outer n surfaces of a circular cylinder Координаты точек x, j (м, рад) Напряжения, МПа, прогиб, см Вариант интерполяционной процедуры Скалярная Векторная Число шагов нагружения 20 30 50 20 30 50 12,77 12,77 12,78 12,67 12,67 12,68 52,82 52,84 52,85 52,53 52,55 52,57 -7,81 -7,81 -7,81 -7,77 -7,77 -7,77 -12,39 -12,40 -12,40 -12,30 -12,31 -12,31 -56,52 -56,54 -56,55 -56,26 -56,27 -56,29 -7,81 -7,81 -7,81 -7,77 -7,77 -7,77 Окончание табл. 2 Координаты точек x, j (м, рад) Напряжения, МПа, прогиб, см Вариант интерполяционной процедуры Скалярная Векторная Число шагов нагружения 20 30 50 20 30 50 12,76 12,77 12,77 12,67 12,68 12,68 52,81 52,83 52,84 52,53 52,55 52,57 -7,81 -7,82 -7,81 -7,77 -7,77 -7,77 -12,40 -12,40 -12,40 -12,30 -12,31 -12,31 -56,53 -56,55 -56,56 -56,26 -56,27 -56,29 -7,81 -7,81 -7,81 -7,76 -7,77 -7,77 w -0,5845 -0,5848 -0,5850 -0,5860 -0,5862 -0,5864 62,12 62,14 62,16 61,74 61,77 61,79 88,38 88,42 88,45 87,11 87,15 87,18 -78,88 -78,91 -78,94 -78,98 -79,01 -79,04 -114,58 -114,63 -114,67 -115,00 -115,05 -115,09 62,11 62,14 62,16 61,74 61,77 61,79 88,38 88,41 88,44 87,12 87,15 87,18 -78,88 -78,92 -78,95 -78,98 -79,02 -79,04 -114,59 -114,64 -114,68 -115,0 -115,05 -115,09 Таблица 3 / Table 3 Численные значения нормальных и сдвиговых напряжений на внутренней и наружней поверхностях эллиптического цилиндра при соотношении a/b=2/1 Numerical values of normal and shear stresses on the inner in and outer n surfaces of theelliptical cylinder with the ratio a/b=2/1 Координаты точек x, j (м, рад) Напряжения, МПа, прогиб, см Вариант интерполяционной процедуры Скалярная Векторная Сетка узлов 21´21 24´24 26´26 28´28 21´21 24´24 26´26 28´28 10,56 10,61 10,55 10,49 8,83 8,60 8,50 8,41 45,49 46,19 46,30 46,36 45,08 44,68 44,52 44,40 -5,54 -6,35 -6,85 -7,34 -4,67 -5,19 -5,55 -5,92 -10,36 -10,42 -10,37 -10,32 -8,63 -8,43 -8,34 -8,26 -47,69 -48,62 -48,87 -49,07 -44,00 -43,85 -43,81 -43,79 -5,54 -6,35 -6,85 -7,33 -4,67 -5,19 -5,55 -5,92 5,29 5,26 5,21 5,17 8,83 8,60 8,50 8,41 26,09 26,44 26,52 26,56 45,08 44,68 44,52 44,40 -2,81 -3,17 -3,40 -3,63 -4,67 -5,19 -5,55 -5,92 -5,14 -5,11 -5,07 -5,03 -8,63 -8,43 -8,34 -8,26 -27,84 28,38 -28,56 -28,70 -44,00 -43,85 -43,81 -43,79 -2,81 -3,17 -3,40 -3,63 -4,67 -5,19 -5,55 -5,92 -0,3333 -0,3433 -0,3464 -0,3487 -0,5551 -0,5579 -0,5592 -0,5602 Таблица 4 / Table 4 Численные значения нормальных и сдвиговых напряжений на внутренней и наружней поверхностях эллиптического цилиндра Numerical values of normal and shear stresses on the inner in and outer n surfaces of the elliptical cylinder Координаты точек x, j (м, рад) Напряжения, МПа, прогиб, см Вариант интерполяционной процедуры Скалярная Векторная Число шагов нагружения 20 30 50 20 30 50 10,48 10,48 10,48 8,49 8,48 8,47 45,77 45,78 45,79 44,19 44,19 44,19 -6,50 -6,50 -6,50 -5,31 -5,30 -5,29 -10,30 -10,30 -10,30 -8,32 -8,31 -8,30 -48,26 -48,27 -48,28 -43,45 -43,45 -43,44 -6,50 -6,50 -6,50 -5,30 -5,30 -5,29 5,19 5,19 5,19 8,49 8,47 8,47 26,25 26,25 26,25 44,19 44,19 44,19 -3,24 -3,24 -3,24 -5,31 -5,30 -5,29 -5,05 -5,05 -5,05 -8,32 -8,31 -8,30 -28,24 -28,24 -28,25 -43,45 -43,45 -43,44 -3,24 -3,24 -3,24 -5,30 -5,30 -5,29 w -0,3498 -0,3499 -0,3501 -0,5671 -0,5674 -0,5677 76,85 76,91 76,95 95,42 95,52 95,61 93,65 93,71 93,75 106,98 107,08 107,15 -91,82 -91,89 -91,95 -106,86 -106,98 -107,08 -120,51 -120,60 -120,67 -130,53 -130,67 -130,78 62,46 62,49 62,51 95,42 95,53 95,61 67,40 67,43 67,46 106,98 107,08 107,16 -76,00 -76,05 -76,08 -106,86 -106,98 -107,08 -91,53 -91,58 -91,62 -130,53 -130,67 -130,78 Исходя из анализа табличного материала, можно сделать вывод, что при определении НДС эллиптического цилиндра в геометрически нелинейной постановке стандартная интерполяционная процедура не приводит к достижению корректного результата. При расчёте эллиптического цилиндра следует применять разработанную векторную форму интерполяции компонент вектора шагового перемещения и шагового вектора углов поворота нормали. Следует отметить, что указанная некорректность стандартной интерполяционной процедуры будет возрастать с увеличением кривизны поперечного сечения эллиптического цилиндра. Выполним расчёт эллиптического цилиндра с вышеприведенными исходными данными при соотношении параметров эллипса поперечного сечения, равного, а также Для наглядности результаты расчетов представим в графической форме. Введём коэффициенты, показывающие отношения нормальных и сдвиговых напряжений в точках и при количестве шагов нагружения, равном 50, и сетке дискретизации 25´25: . На рис. 2-7 показаны изменения вышеупомянутых коэффициентов и в зависимости от соотношения параметров эллипса поперечного сечения эллиптического цилиндра а / b в первом и во втором вариантах расчёта соответственно. Как видно из рисунков, с увеличением отношения а / b значения коэффициентов , в первом варианте существенно возрастают, что является совершенно неприемлемым. Во втором варианте наблюдаем стабильное значение коэффициентов и , равное единице, что и следует ожидать, исходя из условия симметрии расчётной схемы эллиптического цилиндра. Для более полной картины НДС эллиптического цилиндра (при соотношении параметров эллипса поперечного сечения ) на рис. 8 и 9 представлены графики изменения нормальных напряжений Рис. 2. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при стандартной интерполяции искомых величин Fig. 2. Change in the values of the coefficients depending on the ratio of the semiaxes of the ellipse а / b of the cross section of the elliptical cylinder with standard interpolation of the sought values Рис. 3. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при стандартной интерполяции искомых величин Fig. 3. Change in the values of the coefficients depending on the ratio of the semi-axes of the ellipse а / b of the cross section of the elliptical cylinder with standard interpolation of the sought quantities вдоль образующих цилиндра (при и ) при сетке узлов 25´25 и количестве шагов нагружения, равном 30. Анализ графического материала показывает, что в первом варианте расчета, т.е. при скалярной интерполяции искомых величин на шаге нагружения, существенные различия между значениями напряжений сохраняются на протяжении всей длины образующих при и . И, напротив, во втором варианте расчета полное совпадение параметров НДС наблюдается вдоль всей длины образующих эллиптического цилиндра при и , что является необходимым Рис. 4. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при стандартной интерполяции искомых величин Fig. 4. Change in the values of the coefficients depending on the ratio of the semiaxes of the ellipse а / b of the cross section of the elliptical cylinder with standard interpolation of the sought quantities Рис. 5. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при векторной интерполяции искомых величин Fig. 5. Change in the values of the coefficients depending on the ratio of the semi-axes of the ellipse а / b of the cross-section of the elliptical cylinder with vector interpolation of the sought quantities критерием корректности разработанного алгоритма при выбранной расчетной схеме. На рис. 10 и 11 представлены графики изменения нормальных напряжений в зависимости от угловой координаты в концевом и срединном сечениях эллиптического цилиндра при различных вариантах интерполяционной процедуры искомых величин на шаге нагружения. Анализ графического материала, представленного на рис. 10, 11, показывает, что распределение напряжений существенно различается по вариантам расчета, особенно в срединной области дуги Рис. 6. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при векторной интерполяции искомых величин Fig. 6. Change in the values of the coefficients depending on the ratio of the semi-axes of the ellipse а / b of the cross section of the elliptical cylinder with vector interpolation of the sought quantities Рис. 7. Изменение значений коэффициентов в зависимости от отношения полуосей эллипса а / b поперечного сечения эллиптического цилиндра при векторной интерполяции искомых величин Fig. 7. Change in the values of the coefficients depending on the ratio of the semi-axes of the ellipse а / b of the cross-section of the elliptical cylinder with vector interpolation of the sought quantities Рис. 8. Графики изменения нормальных напряжений вдоль образующих цилиндра ( и ) при различных вариантах интерполяции искомых величин Fig. 8. Graphs of changes in normal stresses along the generatrix of the cylinder ( и ) with different options for interpolation of the sought values Рис. 9. Графики изменения нормальных напряжений вдоль образующих цилиндра ( и ) при различных вариантах интерполяции искомых величин Fig. 9. Graphs of changes in normal stresses along the generatrix of the cylinder ( и ) with different options for interpolation of the sought values Рис. 10. Эпюры нормальных напряжений в срединном сечении эллиптического цилиндра при различных вариантах интерполяции искомых величин Fig. 10. Diagrams of normal stresses in the middle section of an elliptical cylinderwith different options for interpolation of the sought values Рис. 11. Эпюры нормальных напряжений в концевом сечении эллиптического цилиндра при различных вариантах интерполяции искомых величин Fig. 11. Diagrams of normal stresses in the end section of an elliptical cylinderwith different options for interpolation of the sought values поперечного сечения и в зоне от до . Так, в первом варианте расчета (при использовании стандартной скалярной интерполяции шаговых искомых величин) в интервале от до наблюдаются скачки в эпюрах напряжений, что явно свидетельствует о некорректности скалярной интерполяционной процедуры, так как в данной зоне отсутствуют локально приложенные нагрузки. Кроме того, на графиках отчетливо прослеживается существенное различие в значениях напряжений в точках и И, напротив, во втором варианте расчета (в котором реализована разработанная векторная форма интерполяции шаговых искомых неизвестных) наблюдается полное совпадение напряжений в точках с угловыми координатами и, что ожидаемо следует из симметрии выбранной расчетной схемы. Анализируя графики напряжений второго варианта расчета, также можно констатировать отсутствие скачков в эпюрах напряжений и симметрию численных значений напряжений, вычисленных на внутренней и наружной поверхностях эллиптического цилиндра относительно нейтральной оси. Заключение На основании вышеизложенного табличного и графического материалов можно сделать вывод о том, что при расчёте эллиптических цилиндров с учетом сдвиговых деформаций в геометрически нелинейной постановке необходимо использовать векторную форму интерполяционной процедуры, применяемую как к компонентам вектора шагового перемещения, так и к компонентам шагового вектора углов поворота нормали (29). Использование стандартной интерполяционной процедуры в задачах исследования НДС эллиптических цилиндров в геометрически нелинейной постановке приводит к получению некорректных результатов.

Об авторах

Ю. В Клочков

Волгоградский государственный аграрный университет

А. Ш Джабраилов

Волгоградский государственный аграрный университет

Т. Р Ищанов

Волгоградский государственный аграрный университет

С. С Марченко

Волгоградский государственный аграрный университет

А. С Андреев

Волгоградский государственный аграрный университет

М. Ю Клочков

Московский государственный университет имени М.В. Ломоносова

Список литературы

Статистика

Просмотры

Аннотация - 270

PDF (Russian) - 229