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

Обложка


Цитировать

Полный текст

Аннотация

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

Полный текст

  1. Предварительные сведения

    Целью исследования является развитие непрямого метода граничных элементов (НМГЭ) или, как его еще называют, метода компенсирующих нагрузок для решения задач деформирования изотропных пластин на упругом основании со сложным контуром в условиях термомеханического нагружения (сосредоточенные, распределенные нагрузки и температурные поля) при различных граничных условиях. Рассматриваются малые деформации тонкой линейно-упругой пластинки на упругом основании, деформирование которой описывается моделью, основанной на гипотезах Кирхгофа — Лява в рамках

    теории малого изгиба [1].

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

    Для решения поставленной задачи МГЭ необходимо предварительно определить фундаментальное решение задачи изгиба изотропной пластины на упругом основании Винклера и Пастернака — Власова. Для определения фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера и Пастернака — Власова и затем с помощью НМГЭ решения этих задач была использована универсальная методика, предложенная в [1–4]. Для отработки методики были решены несколько задач изгиба изотропных пластин на упругом основании Винклера и Пастернака — Власова.

     

  2. Постановка задачи

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

    αt

    D22W = pz (x, y) p(W ) (1 + ν)

    T,

    h

    где 2 — оператор Лапласа; W — прогиб точки срединной поверхности пластины; pz — интенсивность

    нормального давления, действующего на пластину; p — реактивное давление, которое зависит от

    упругого основания, на котором лежит пластина; ν — коэффициент Пуассона (коэффициент поперечной деформации); αt — коэффициент линейного температурного расширения; h — толщина пластины; T — изменение температуры между поверхностями пластины.

    Рассмотрим упругое основание, которое характеризуется следующей функцией:

    p(W ) = K1W + K2W 2 + K3W 3 g1

    2W

     

    image

    ∂x2

    g2

    2W

    image

    ∂y2 ,

    где Ki (i = 1, 2, 3) и gi (i = 1, 2) — коэффициенты упругого основания. В таблице приведено несколько возможных вариантов оснований.

     

     

    Наиболее распространенные варианты упругих оснований The most common options for elastic foundations

    Таблица

    Table

     

    K1 ̸= 0 ; K2 = K3 = g1 = g2 = 0

    Линейное изотропное основание Винклера (Циммермана — Фусса)

    K1 ̸= 0; K2 = K3 = 0; g1 = g2 = g ̸= 0

    Линейное основание Пастернака — Власова

    K1 ̸= 0; K2 = K3 = 0; g1 ̸= g2 ̸= 0

    Линейное ортотропное основание

    K1 ̸= 0; K2 ̸= 0; K3 ̸= 0

    Нелинейное основание

     

    В дальнейшем будем рассматривать линейные упругие основания Винклера и Пастернака — Власова как самые распространенные.

    Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

    Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 35

     

  3. Определение фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера и Пастернака — Власова

    Пусть дана система линейных дифференциальных уравнений в частных производных [5]:

    L0U (x, y) = F (x, y), (1)

    image

    image

    где L0 = [Llm( , )] — исходный линейный дифференциальный оператор в частных производных,

    ∂x ∂y

    image

    image

    l, m = 1, N ; U (x, y) = [ul(x, y)]T — векторная функция, подлежащая определению; F (x, y) = [fl(x, y)]T — заданная векторная функция правых частей; l, m = 1, N .

    Решение системы (1) U (x, y) представимо в виде свертки:

    0

     

    U (x, y) = L1F (x, y) = G F =

    G(x ξ, y η)F (ξ, η) dξdη, ξ, η ,

    0

     

    где L1

    • интегральный оператор, ядром которого является матрица фундаментального решения

      G(x ξ, y η) системы линейных дифференциальных уравнений; — область определения

      дифференциального оператора L0 .

      Матрица фундаментального решения G(x ξ, y η) определяется из выражения вида:

      L0G(x ξ, y η) = δ(x ξ, y η)I, (2) где δ(x ξ, y η) — двумерная дельта-функция Дирака; I — единичная матрица размерностью N × N . Фундаментальные решения определяются с точностью до решения однородной системы уравнений

      L0G(x ξ, y η) = 0 и, кроме того, являются обобщенными функциями.

      Из (2) видно, что фундаментальное решение G(x ξ, y η) зависит только от свойств

      дифференциального оператора L0.

      В дальнейшем фундаментальное решение G(xξ, y η) будем обозначать G(t, ζ) = G(tζ), где t(x, y)

      и ζ(ξ, η).

      Элементы матрицы G(t, ζ) определятся из решения систем дифференциальных уравнений вида:

       

      где δlm — символ Кронекера.

      N

      image

      Lmk Gkl(t, ζ) = δlmδ(t ζ), (l, m = 1, N ), (3)

      k=1

      Решение систем (3) осуществляется, например, методом двумерного интегрального преобразования

      Фурье [1–4]. Без ограничения общности и в целях экономии места далее предполагаем, что точка приложения единичных сосредоточенных нагрузок, которые моделируются обобщенной δ-функцией Дирака, находится в начале координат.

      image

      Двумерная трансформанта (образ, изображение) f (ξ, η) функции f (x, y) определяется выражением вида:

       

      ++

      1

      image

      F [f (x, y)] = f (ξ, η) =

      2π

      image

      f (x, y)ei(ξx+ηy)dxdy, i = 1.

      −∞ −∞

      image

      Соотношение между трансформантами производных от f (x, y) и самой трансформантой f (ξ, η) может быть определено следующим выражением:

      (n+m)f (x, y)

      image

      image

      F [ ∂xn∂ym

      ](ξ, η) = ()n()m f (ξ, η).

      Трансформанта двумерной обобщенной δ-функции Дирака, которая, в силу присущих ей свойств, представима в виде δ(x, y) = δ(x)δ(y), определяется из следующего соотношения: F [δ(x, y)] =

      ++

      =

       

      1 ∫ ∫

      image

      2π

      image

      2π

       

      δ(x, y)ei(ξx+ηy)dxdy = 1 .

      −∞ −∞

      Трансформанты для производных δ-функции Дирака примут вид:

      image

      − −

       

      F [δ(m+n)(x, y)] = 1 ( )m( )n.

      2π

      Системы (3) в пространстве трансформант примут вид:

      N

      image

      Lmk (iξ, )Gkl(ξ, η) =

      k=1

      image

      1

      image

      2π δlm, (l, m = 1, N ).

      Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

      36 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

       

      Формула обращения двумерного интегрального преобразования Фурье (прообраз, оригинал) имеет

      ++

      image

      2π

       

      вид: f (x, y) = 1

      f (ξ, η)ei(ξx+ηy)dξdη.

      −∞ −∞

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

      ++

      (n+m)f (x, y)

      image

      1 ∫ ∫ n m

      i(ξx+ηy)

      image

      image

      ∂xn∂ym = 2π

       

      −∞ −∞

      ()

      ()

      f (ξ, η)e

      dξdη.

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

      Q(z) — функция, где z — комплексное число, удовлетворяющая условиям:

      1. аналитическая при Im z ): 0 всюду, за исключением числа полюсов, лежащих на действительной оси;

      2. не имеет полюсов на действительной оси;

      3. при z → ∞ функция zQ(z) 0 равномерно при 0 arg z π;

      4. интегралы

       

      Q(x)dx и

      0

      0

      Q(x)dx являются сходящимися, то

      −∞

      +

       

      Кроме того,

      Q(x)dx = 2πi

      −∞ Imz�0

       

      Res Q(z).

      Q(x) cos mx dx = πi

      0 Imz�0

      Res Q(z) eim z ;

       

      Q(x) sin mx dx = π

      0 Imz�0

       

      Res Q(z) eim z.

      Если z = a — полюс функции Q(z) порядка m, то вычет функции Q(z) определяется по формуле:

      1

      image

      Res Q(a) =

       

      lim

      image

      dm1

      [(z a)m Q(z)].

      (m 1)! za dzm1

      Для вычисления расходящихся интегралов далее будем вводить в рассмотрение обобщенную функцию Θ(η), производная от которой также будет обобщенной, и пользоваться правилом дифференцирования обобщенных функций [4].

      Рассмотрим вычисление несобственных интегралов следующих типов:

      ++

      1

      image

      I1(x, y) = 2π

      1

      image

      (ξ2 + η2) e

      i(ξx+ηy)

      dξdη,

      −∞ −∞

       

      ++

      1

      image

      I2(x, y) = 2π

      1

      image

      (ξ2 + η2 + λ2) e

      i(ξx+ηy)dξdη.

      −∞ −∞

      Интеграл I1(x, y) определим, воспользовавшись свойством интегралов от четных и нечетных функций на интервале, симметричном относительно начала координат:

      ++

      1

      1

      i(ξx+ηy) 2

      cos ξxdξ

      1 ηx

      I1(x, y) = 2π

      −∞ −∞

      (ξ2 + η2) e

      dξdη =

      π

      0

      cos ηydη

      0

      (ξ2 + η2) = η e

      0

      cos ηydη,

      где была использована теория вычетов:

      cos ξx

      eizx

      eiz x 1 π 1

      image

      image

      1

       

      1 1 ηx

      (ξ2 + η2)= πi

      Res = πi

      (z2 + η2) 2z

      = e.

      2 η

      0 Imz�0

      1 1z1 =

      Внешний интеграл вычислим с помощью теории обобщенных функций [4]. Рассмотрим обобщенную функцию вида:

      Θ(η) = ln η+ =

      { 0, η 0,

      ln η, η > 0.

      image

      η

       

      Производная этой функции есть обобщенная функция Θ(η) = 1 , η > 0.

      Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

      Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 37

       

      Согласно правилу дифференцирования обобщенных функций, имеем

      , φ) = , φ) =

      φ(η)Θ(η)d(η),

       

      где φ(η) – основная функция (бесконечное число раз дифференцируемая и в бесконечно удаленной точке сама функция и ее бесконечное число производных стремятся к нулю).

      Здесь в качестве основной функции выступает функция φ(η) = eηx cos ηy. Учитывая все вышесказанное, получим:

       

      I1(x, y) =

      1

      eηx

      η

      0

      cos ηydη = x

      0

      1 2

       

      eηx

       

      2

       

      cos ηy ln η dη + y

       

      eηx

      0

       

      sin ηy ln ηdη =

      image

      = 2 ln(x

      + y ) C = (ln r + C),

      где входящие интегралы являются табличными; r2 = x2 + y2, C = 0.5772157 ... – постоянная Эйлера.

      Рассмотрим интеграл I2(x, y):

      ++

      1

      image

      I2(x, y) = 2π

      1

      image

      (ξ2 + η2 + λ2) e

      image

      i(ξx+ηy)dξdη = 2

      π

      cos ηydη

      cos ξxdξ

      image

      (ξ2 + η2 + λ2) .

      −∞ −∞ 0 0

      Аналогично рассмотренному выше

      cos ξx π 1

      отсюда

      image

      xη2 +λ2

      image

      (ξ2 + η2 + λ2) = 2 η2 + λ2 e ,

      0

       

      I2(x, y) =

      image

       

      1 xη2 +λ2

      η2 + λ2 e

      0

       

      cos ηy dη = K0(λr),

      где K0(λr) — модифицированная функция Бесселя второго рода нулевого порядка; этот интеграл также является табличным.

      Опираясь на все вышеизложенное, можно показать, как получаются фундаментальные решения задачи изгиба изотропной пластины на упругом основании Винклера и Пастернака — Власова.

       

      3.1. Нахождение фундаментального решения с помощью метода последовательного интегрирования

      Запишем, например, уравнение изгиба изотропной пластины на упругом основании Винклера в виде:

      4 4

      ( + 2

      4

       

      + + χ4)W (x, y) = q(x, y)

      q(x, y)

      ( L W (x, y) = ),

      image

      ∂x4

      image

      ∂x2∂y2

      image

      image

      ∂y4

      image

      D 0 D

      image

      image

      3

       

      где χ4 = K1 ; D = Eh

    • изгибная (цилиндрическая) жесткость пластины.

    D 12(1ν2 )

    Представим дифференциальный оператор изгиба пластины на упругом основании L0 в виде двух

    последовательных операторов Гельмгольца:

    4 4 4 4

    2 2

    2 2 2 2

    image

    image

    image

    image

    image

    image

    image

    L0 = ( ∂x4 + 2 ∂x2y2 + y4 + χ ) = ( ∂x2 + y2 )( ∂x2 + y2 + ).

    Фундаментальное решение ищем в виде:

     

    L0G1(x, y) =

    δ(x, y)

    image

    .

    D

    После введения тривиального обозначения получим

    2

    image

    ( x2 +

    2

     

    image

    ∂y2

    2) G˜1(x, y) =

    δ(x, y)

    image

    .

    D

    Далее к этому дифференциальному уравнению последовательно применяем двумерное интегральное преобразование Фурье, выделение трансформанты G˜1(x, y), применение формулы обращения двумерного интегрального преобразования Фурье. В результате G˜1(x, y) примет вид:

    image

    1

     

    G˜1(x, y) = 2πD K0(

    image

    iχr).

    Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

    38 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

     

    Далее по известной величине G˜1(x, y) восстанавливаем G1(x, y):

    i G1(x, y) = 2D

    0

    image

    image

    rK0(

    0

     

    (r0)H(2)(

    (r r0))dr0,

    0

     

    где H(2)

    — функция Ханкеля 0-го порядка 2-го рода.

    Выполняя аналогичные действия с дифференциальным оператором изгиба пластины на упругом основании L0 в виде двух переставленных последовательных операторов Гельмгольца, получим

    4 4 4 4

    2 2

    2 2 2 2

    image

    image

    image

    image

    image

    image

    image

    L0 = ( ∂x4 + 2 ∂x2y2 + y4 + χ ) = ( ∂x2 + y2 + )( ∂x2 + y2 ),

    2

    image

    ( ∂x2

    2

    image

    + ∂y2

    + 2)G˜2(x, y) =

    i

    δ(x, y)

    image

    ,

    D

    image

    G˜2(x, y) =

    image

    0

     

    H(2)(iχr),

    4D

    i G2(x, y) = 2D

    0

    image

    0

     

    rH(2)(

    (r0)K0(

    image

    (r r0))dr0.

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

    l2 r

     

    image

    где χ4 = 1

    image

    image

    G(x, y) = G1(x, y) + G2(x, y) = 2πD kei0( l ),

    image

    image

    = K1 ; kei0( r ) — функция Томпсона-Кельвина.

    l4 D l

    Для подтверждения правильности найденного фундаментального решения задачи изгиба изотропной

    пластины на упругом основании Винклера рассмотрим еще 2 варианта его нахождения:

    а) трансформанта Фурье фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера имеет вид:

    1

    image

    G(ξ, η) =

    1 i 1

    image

    image

    = (

    1

    image

    ).

    2πD (ξ2 + η2)2 + χ4

    4πDχ2

    ξ2 + η2 + 2

    ξ2 + η2 2

    Восстанавливая его оригинал с точностью до решения однородного уравнения, получим

    G(x, y) =

    i

    image

    4πDχ2

    image

    K0(iχr)

    1

    image

    82

    image

    0

     

    H(2)(iχr) =

    l2

    image

    2πD

    image

    r kei0( l );

    б) трансформанта Фурье фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера имеет вид:

    1

    image

    G(ξ, η) =

    1 1 1

    image

    image

    = Im .

    2πD (ξ2 + η2)2 + χ4

    Восстанавливая оригинал, получим:

    2πDχ2

    ξ2 + η2 + 2

    1

    image

    l2 r

    G(x, y) = 2πDχ2 Im(K0(χr

    image

    i)) = 2πD kei0( l ).

    Таким образом, различные варианты получения фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера совпали.

    Для проверки правильности полученного фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера были использованы две методики: с помощью формулы дифференцирования обобщенных функций [4]

    Φ

    image

    ∂xi

    Φ

    image

    = (

    ∂xi

    )G + (1)i1δ(x1, . . . , xn)

    Γ

     

    Φdx1 . . . dxi

     

    1dx

     

    i+1

     

    . . . dxn,

    image

    ∂xi

     

    где ( Φ )G — обычная производная от функции Φ, Γ — граница области G (один из контуров, внутри которого находится особенность) и с помощью проверки равновесия пластины на упругом основании Винклера, ограниченной кривой, при действии на нее единичной нагрузки. Результаты проверок подтверждают правильность найденного фундаментального решения задачи изгиба изотропной пластины на упругом основании Винклера.

    Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

    Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 39

     

    Все вышеприведенные действия были также проведены и для получения фундаментального решения задачи изгиба изотропной пластины на упругом основании Пастернака — Власова, изгиб которой описывается дифференциальным уравнением вида:

     

    L0W (x, y) =

    pz (x, y) , D

    где L0 = 2 2p2∆+χ4 — дифференциальный оператор; kz, kt — параметры упругого основания (первый и второй коэффициенты постели или коэффициенты сжатия и сдвига соответственно);

    2

    image

    2p2 = g

    kth

    image

    =

    image

    image

    image

    , χ4 = K1 = kz = 1 .

    D 4D

    D D l4

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

    image

    4πs2 D

     

    других, которые получаются из них как частные случаи: a) p > χ G(x, y) = 1

    image

    [K0(p2 + s2r)

    image

    image

    K0(p2 s2r)], s2 = p4 χ4;

    2

    image

    image

    b) p < χ G(x, y) = 1

    image

    image

    kei0(r/l, φ), β2 = χ4 p4, φ = 1 arctanβ ;

    2πβ2 D

    image

    4πDχ2

     

    c) p = χ G(x, y) = 1

    image

    4πp2 D

     

    d) χ = 0 G(x, y) = 1

     

    0

     

    kei (r/l, 0);

    image

    [ln r + K0(2pr)];

    2 p2

    image

    2

     

    e) p = 0 G(x, y) = l

    kei0(r/l);

    2πD

    image

    8πD

     

    f) p = χ = 0 G(x, y) = 1

     

    r2ln r.

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

    Для реально существующих материалов пластин, лежащих на упругом основании Пастернака — Власова, имеет место случай b):

    1

    image

    G(x, y) = 2πβ2D kei0(r/l, φ),

    где kei0(x, φ) — обобщенная функция Томпсона-Кельвина (мнимая часть обобщенной модифицированной функции Бесселя второго рода нулевого порядка K0).

     

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

    Рассмотрим тонкую линейно-упругую пластину, срединная поверхность которой занимает область +, ограниченную гладким кусочно-ляпуновским контуром Γ (рис. 1). Дополним область + областью до бесконечной области. По контуру Γ к бесконечной пластине приложим компенсирующие нагрузки q(ζ), m(ζ). Нагрузки q(ζ) и m(ζ) – распределенные по контуру Γ усилие, нормальное к поверхности пластины, и момент вокруг касательной к контуру Γ соответственно.

     

    image

     

    Рис. 1. Тонкая линейно-упругая пластина Fig. 1. Thin linear elastic plate

    Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

    40 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

     

    Решение задачи изгиба многосвязной пластины рассматривается как сумма основного (частного) и компенсирующего решений [6]:

    M

     

    W (t) = Wr (t) + [

    i=1 Γi

    G(t, ζ)qi(ζ)

    i

     

    ∂G(t, ζ)

    image

    m (ζ)ds(ζ)], (4)

    ∂n1

    где Wr (t) =

    +

    n

    G(t, ζ)pz (ζ)dΩ(ζ) + G(t, ζi)Pi(ζi).

    i=1

    Здесь t(x, y) — точка внутренней области +; ζ(ξ, η)) — точка контура Γ; G(t, ζ) — фундаментальные

    решения задачи изгиба пластины, лежащей на упругом основании Винклера или Пастернака — Власова; pz (ζ) – интенсивность нормального давления, действующего на пластину; Pi — модуль сосредоточенной силы, приложенной в точке ζi(i = 1, 2, ..., n).

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

    На контуре пластины рассмотрим следующие граничные условия:

    ∂W

     

    W = 0, W,n =

    image

    = 0 жесткая заделка, (5)

    ∂n

    W = 0, Mn = 0 шарнирное опирание, (6)

    Vn = k1W, Mn = k2W,n упругое закрепление. (7)

    Здесь Mn, Vn — изгибающий момент и обобщенная поперечная сила на контуре Γ; k1, k2 - постоянные, определяющие упругие свойства закрепления (при k1 = k2 = 0 получим граничные условия свободного края).

    Усилия и моменты на контуре пластины в глобальной декартовой системе координат X, Y

    определяются по формулам [4]:

     

    M =

    Mn = Mx cos2 α + My sin2 α + Mxy sin 2α;

    MT = Mx sin2 α + My cos2 α Mxy sin 2α;

    Mx My sin 2α + M cos 2α; Q = Q cos α + Q sin α;

    2 xy n x y

    ∂M

    Vn = Qn + s ,

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

    Mx = D(W, xx +νW, yy ); My = D(W, yy +νW, xx ); Mxy = D(1 ν)W, xy ;

    ∂ ∂

    image

    image

    Qx = D ∂x W = D(W, xxx +W, xyy ); Qy = D y W = D(W, yyy +W, xxy ).

    Для изотропной пластины, лежащей на упругом основании, моменты и усилия на контуре пластины можно также определить в локальной системе координат n, τ , построенной в точке t(x, y) контура, где записываются граничные условия. При этом ось n — внешняя нормаль к контуру в этой точке; а ось τ — касательная к контуру в этой точке.

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

    Mn = D(∆W (1 ν)W, ττ ); Mτ = D(∆W (1 ν)W, nn );

    ∂M

    M = D(1 ν))W, ; Qn = DW,n Vn = Qn + s ,

    2 2

    image

    ∂n2

     

    +

     

    где ∆ = 2 =

    image

    ∂τ 2

    • оператор Лапласа в локальной системе координат.

      image

      ∂s

       

      При определении производной по дуговой координате контура ∂M

      следует учитывать, что α =

      image

      image

      ∂s

       

      = ±

       

      ρ

       

      = α(s); k = α 1

    • кривизна контура в рассматриваемой точке (алгебраическая величина), ρ

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

    ∂ ∂3W

    image

    image

    Vn = D( n W + (1 ν)[ ∂n∂τ 2 k(2W,xy sin 2α + [W,xx W,xy ] cos 2α)]) =

    Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

    Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 41

     

    ∂ ∂3W

    2W

    2W

    image

    image

    image

    = D( n W + (1 ν)[ ∂n∂τ 2 k( n2

    image

    ∂τ 2 )]).

    Компенсирующие нагрузки определяются из граничных условий (5)–(7) на контуре пластины, которые записываются на контуре Γ для области +.

    Известно, что выражение (4) в точности удовлетворяет дифференциальным уравнениям задачи изгиба изотропной пластины на упругом основании Винклера и Пастернака — Власова. Чтобы

    (4) было решением краевой задачи, необходимо из системы сингулярных интегральных уравнений, которая получается при подстановке (4) в краевые условия (5)–(7) на контуре пластины, определить компенсирующие нагрузки q(ζ), m(ζ). При этом нужно совершить предельный переход точки t(x, y) из внутренней области + на границу Γ. Тогда ядра интегральных уравнений будут определены во всех точках контура, за исключением точки, где t(x, y) = ζ(ξ, η). В этой точке ядра потенциалов будут иметь особенности.

    В области + изгибающие и крутящий моменты, а также поперечные силы для изотропной пластины, лежащей на упругом основании, по известным компенсирующим нагрузкам q(ζ), m(ζ), определяются по формулам:

    2G

    2G

    3G

    3G r

    Mx(t) = D

    Γ

    image

    image

    image

     

    [( + ν )q(ζ) (

    ∂x2 ∂y2 ∂x2∂n1

    image

    + ν

    ∂y2∂n1

    )m(ζ)]ds(ζ) + Mx (t);

    2G

    2G

    3G

    3G r

    My (t) = D

    Γ

    image

    image

    image

     

    [( + ν )q(ζ) (

    ∂y2 ∂x2 ∂y2∂n1

    image

    + ν

    ∂x2∂n1

    )m(ζ)]ds(ζ) + My (t);

    2G

    3G r

    Mxy (t) = D(1 ν)

    Γ

    image

    image

     

    [ q(ζ)

    ∂x∂y ∂x∂y∂n1

    m(ζ)]ds(ζ) + Mxy (t);

    G

    2G r

    Qx(t) = D [

    Γ

    image

     

    q(ζ)

    ∂x ∂x∂n1

    m(ζ)]ds(ζ) + Qx(t);

    G

    2G r

     

    где

    Qy (t) = D [

    Γ

    image

     

    q(ζ)

    ∂y ∂y∂n1

    m(ζ)]ds(ζ) + Qy (t),

    2 r 2 r

    2 r 2 r 2 r

    ∂ W (t)

    Mr

    ∂ W (t) r

    ∂ W (t)

    ∂ W (t) r

    ∂ W (t)

    x (t) = D(

    x2 + ν

    y2 ); My (t) = D(

    y2 + ν

    ∂x2 ); Mxy (t) = D(1 ν)

    ;

    ∂x∂y

    3 r 3 r

    3 r 3 r

    ∂ W (t)

    Qr

    ∂ W (t) r

    ∂ W (t)

    ∂ W (t)

    x(t) = D(

    ∂x3 +

    ∂x∂y2 ); Qy (t) = D(

    ∂x2∂y +

    y3 ).

     

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

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

    С учетом предельных значений потенциалов по методике [7], а также результатов дифференцирования фундаментального решения задачи изгиба пластины на упругом основании запишем сингулярные интегральные уравнения, из решения которых определяются компенсирующие нагрузки q(ζ), m(ζ) для различных граничных условий. При этом предельные значения потенциалов берутся в области +.

    1. Жесткая заделка

      ∫ ∫

      G(t, ζ)q(ζ)ds(ζ)

      Γ Γ

      ∂G(t, ζ)

      ∂n1

      m(ζ)ds(ζ) + W

      r (t) = 0,

      ∂G(t, ζ)

       

      image

      ∂n

      Γ

      q(ζ)ds(ζ)

      Γ

      2G(t, ζ)

       

      image

      ∂n∂n1

       

      m(ζ)ds(ζ) +

      ∂Wr (t)

       

      image

      ∂n

       

      = 0.

      Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

      42 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

       

    2. Шарнирное закрепление

       

      Γ

       

      G(t, ζ)q(ζ)ds(ζ)

      Γ

       

      ∂G(t, ζ)

      ∂n1

       

      m(ζ)ds(ζ) + W

       

      r (t) = 0,

      1

      image

      2 m(t) D

      Γ

      [∆G(t, ζ) (1 ν)

      2G(t, ζ)

       

      image

      ∂τ 2

      ]q(ζ)ds(ζ)+

       

    3. Свободный край

    +D [

    Γ

    G(t, ζ)

     

    image

    ∂n1

    (1 ν)

    3G(t, ζ)

     

    image

    ∂τ 2∂n1

    n

     

    ]m(ζ)ds(ζ) + Mr (t) = 0,

    1

    image

    2 m(t) D

    Γ

    [∆G(t, ζ) (1 ν)

    2G(t, ζ)

     

    image

    ∂τ 2

    ]q(ζ)ds(ζ) + D

    Γ

    G(t, ζ)

    image

    [

    ∂n1

    (1 ν)

    3G(t, ζ)

     

    image

    ∂τ 2∂n1

     

    ]m(ζ)ds(ζ)+

    1

    image

    2 q(t) D

    Γ

     

    G(t, ζ)

    image

    [

    ∂n

     

    + (1 ν)[

     

    3G(t, ζ)

     

    image

    ∂n∂τ 2

    n

     

    +Mr (t) = 0,

    2G(t, ζ)

    image

    k(t)( n2

     

    2G(t, ζ)

     

    image

    ∂τ 2

    )]]q(ζ)ds(ζ) + D

    Γ

     

    2G(t, ζ)

    image

    [ +

    ∂n∂n1

    4G(t, ζ)

    image

    1

     

    +(1 ν)[ ∂n∂τ 2n

    3G(t, ζ)

    image

    1

     

    k(t)( n2n

    3G(t, ζ)

    image

    ∂τ 2∂n1

    n

     

    )]]m(ζ)ds(ζ) + V r (t) = 0,

    где введены обозначения

    Mr r r

    n(t) = D[W,nn (t) + νW,ττ (t)];

    3 r 2 r 2 r

    V r r

     

    ∂ W (t)

    ∂ W (t)

    ∂ W (t)

    n (t) = D( n W

    (t) + (1 ν)[

    ∂n∂τ 2 k(t)(

    ∂n2

    ∂τ 2 )]).

    Функции плотностей q(ζ) и m(ζ) удовлетворяют условию Гельдера:

    |q(t) q(ζ)| C|t ζ|α; |m(t) m(ζ)| C|t ζ|α,

    где C = const > 0, 0 < α 1.

    Внеинтегральные члены q(t)/2, m(t)/2 в интегральных уравнениях получаются на основе анализа предельных значений потенциалов [7].

     

  6. Численная реализация

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

    N

    ΓN =

    Γj и ΓN Γ при N → ∞ или max |Γj | → 0, где N – количество элементов контура ΓN .

    j=1 j

    В пределах элемента компенсирующие нагрузки q(ζ) и m(ζ) считаются постоянными. Точки коллокации берутся в узлах, которые располагаются в серединах элементов Γj .

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

    N N

    image

    q(ζj )Gij m(ζj )Hij + Ai = 0, 1, N,

    j=1

    N

    j=1

    N

    image

    q(ζj )Rij m(ζj )Mij + Bi = 0, 1, N,

    j=1

    j=1

    где коэффициенты Gij, Hij, Rij, Mij в зависимости от вида граничных условий принимают соответствующие выражения и представляют собой интегралы по элементу Γj ; Ai, Bi — правые части системы уравнений, принимающие соответствующие выражения в зависимости от граничных условий; q(ζj ) = qi, m(ζj ) = mj — искомые значения компенсирующих нагрузок, принимающие постоянные значения в пределах каждого элемента Γj . Система состоит из 2N уравнений с 2N неизвестными, т. е. является замкнутой и имеет единственное решение.

    Интегралы Gij, Hij, Rij, Mij при i ̸= j не имеют особенностей и вычисляются по восьмиузловой

    квадратурной формуле Гаусса. При i = j ядра интегралов имеют особенности типа ln r, 1/r, 1/r2

    при r 0 и, кроме того, необходимо учесть внеинтегральные члены уравнений. Такие сингулярные

    интегралы по Γj вычисляются аналитически.

    Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

    Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 43

     

  7. Вычисление сингулярных интегралов по элементам контура

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

    особенности при x ξ, y η. Считаем, что компенсирующие нагрузки в пределах этих элементов

    постоянны, поэтому они входят в интегральные уравнения в виде множителей перед интегралами от

    фундаментального решения и его производных.

     

    image

     

    Рис. 2. Прямолинейный граничный элемент Fig. 2. Straight boundary element

     

    Пусть Γj (j = 1, 2, ..N ) — один из элементов контура Γ, представляющий собой отрезок прямой (рис. 2), соединяющей точки (xA, yA) и (xB, yB ), а t(x, y) — середина отрезка Γj :

    x = (xA + xB )/2; y = (yA + yB )/2.

    Введем координату t0, определяющую положение точки ζ(ξ, η) на Γj :

     

    где t0[1; +1].

    ξ = (xA + xB )

    2

    + (xB xA)

    2

     

    t0; η =

    (yA + yB ) +

    2

    (yB yA)

    2

     

    t0,

    Расстояние между точками t и ζ равно

    image

    image

    image

    2

     

    B

     

    r = (x ξ)2 + (y η)2 = |t0| (x

     

    • xA

       

      )2 + (yB

       

    • yA

    image

    )2 = |t0| δ,

    2

    где δ — длина отрезка Γj .

    Фундаментальное решение и его производные зависят от аргументов x ξ, y η, которые

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

    x ξ =

    0 − −

     

    (xB xA) t ; y η = 2

    (yB yA) t .

    2 0

    Выразим элемент длины ds через координату t0

    image

    ds = ξ )2 + (η )2dt0 =

    image

    (

     

    (xB xA)

     

    )2 + (

    (yB yA)

    δ

    )2dt0 =

     

    dt0.

    ( t0 t0

    2 2 2

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

     

  8. Результаты решения некоторых задач

Пример 1

Кольцевая пластина с внешним радиусом R = 0.1 м и внутренним радиусом r = 0.05 м лежит на упругом основании Пастернака — Власова с коэффициентами постели kz = 15400 MH/м3, kt =

= 10000 МН/м3 (рис. 3). Внутренний контур шарнирно оперт, а внешний — жестко защемлен.

На пластину действует постоянная равномерно распределенная поперечная нагрузка pz = 1 МПа. Физические параметры: h = 0.01 м, E = 2 · 105 МПа, ν = 0.3. При численном решении каждый контур

image

k

 

разбивался на 40 элементов, а частное решение принято в виде Wr (t) = pz . На рис. 4, а, б приведены

z

результаты решения поставленной задачи: графики распределения прогибов W , углов поворота Wr ,

радиального Mr и тангенциального Mφ моментов по радиусу пластины. Точки соответствуют численному решению НМГЭ, а непрерывные линии — точному решению задачи, которое ищется в виде:

r r r r r

W (r) = C1 ber( l , φ) + C2 bei( l , φ) + C3 ker( l , φ) + C4 kei( l , φ) + W

(r),

Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

44 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

 

где постоянные Ci определяются из граничных условий на контурах пластины. Как видно, в пределах точности графиков результаты решений совпадают.

 

image

 

Рис. 3. Кольцевая пластина для расчета Fig. 3. Ring plate for calculation

 

image

 

a б

 

Рис. 4. Результаты расчета Fig. 4. Calculation results

 

Пример 2

 

 

image

 

Рис. 5. Пологая сферическая оболочка Fig. 5. Shallow spherical shell

Вестник Самарского университета. Естественнонаучная серия. 2021. Том 27, № 2. С. 33–47

Vestnik of Samara University. Natural Science Series. 2021, vol. 27, no. 2, pp. 33–47 45

 

Квадратная в плане пологая изотропная сферическая оболочка радиуса R с угловым вырезом размером a с граничными условиями шарнирного закрепления (пунктирная линия), жесткой заделки (наклонная линия) и свободного края (нет линии) (рис. 5) лежит на упругом основании Винклера. Оболочка находится под действием постоянного температурного поля, которое по толщине пластины изменяется по линейному закону, и сосредоточенной силы P , приложенной в начале координат. Перепад

температуры по толщине пластины T = 100 град. Коэффициент температурного расширения материала примем следующим: αt = 0.125 · 104 град1. Геометрические и механические параметры пластины: a = 0.1 м, h = 0.01 м, R = 0.5 м, Kz = 15400 МН/м3, P = 4000 Н, E = 2 · 105 МПа, ν = 0.3.

Контур разбит на 60 элементов, по 10 на каждую сторону.

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

2 2 kth2 2 1 Eh

pz (x, y)

αt 2

image

∇ ∇ W

4D W + D (kz + R2 )W =

D (1 + ν) h T,

конкретизируя его на случай основания Винклера (kt = 0).

На рис. 6 представлены результаты решения поставленной задачи: на рис. 6, а прогиб w, а на рис. 6, б напряжения (как нормальные, так и касательные) σx, σy, τxy при z = h/2 в сечении пластины x = 0. На графиках напряжений сплошная линия соответствует σx; пунктирные линии (линия с большей длиной пунктира и линия с меньшей длиной пунктира) – соответственно σy и τxy .

 

image

 

a б

 

Рис. 6. Результаты расчета Fig. 6. Calculation results

 

Пример 3

 

image

 

a б

 

Рис. 7. Двусвязная пластина сложной формы Fig. 7. Doubly connected plate of complex shape

Великанов П.Г., Куканов Н.И., Халитова Д.М. Использование непрямого метода граничных элементов...

46 Velikanov P.G., Kukanov N.I., Khalitova D.M. Use of the indirect boundary elements method for isotropic plates...

 

Двусвязная пластина сложной формы с различными граничными условиями на отдельных участках лежит на упругом основании и находится под действием термомеханического нагружения (сосредоточенной силы P , приложенной в начале координат и линейного по толщине температурного

поля) (рис. 7, а). Параметры пластины и поля: h = 0.01 м, a = 0.1 м, E = 2 · 105 МПа, ν = 0.3, kz =

= 15400 МН/м3, P = 5000 Н, T = 100 град, αt = 0.125 · 104 град1.

На рис. 7, б приведена поверхность распределения прогиба W с линиями уровня.

 

Выводы

Таким образом, в данной статье сначала различными методами были получены фундаментальные решения, а затем НМГЭ были решены задачи расчета изотропных пластин на упругом основании Винклера и Пастернака — Власова. Для отработки методики были решены тестовые задачи, подтвердившие высокую точность результатов.

×

Об авторах

П. Г. Великанов

Казанский (Приволжский) федеральный университет

Автор, ответственный за переписку.
Email: pvelikanov@mail.ru
ORCID iD: 0000-0003-0845-2880

кандидат физико-математических наук, доцент кафедры теоретической механики

Россия, 420008, Российская Федерация, г. Казань, ул. Кремлевская, 18

Н. И. Куканов

Ульяновский государственный технический университет

Email: kukanov_n_i@mail.ru
ORCID iD: 0000-0003-0880-4591

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

Россия, 432027, Российская Федерация, г. Ульяновск, ул. Северный Венец, 32

Д. М. Халитова

Казанский (Приволжский) федеральный университет

Email: diana982000@gmail.com
ORCID iD: 0000-0002-2239-9222

магистрант

Россия, 420008, Российская Федерация, г. Казань, ул. Кремлевкая, 18

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

  1. Великанов П.Г. Метод граничных интегральных уравнений для решения задач изгиба изотропных пластин, лежащих на сложном двухпараметрическом упругом основании // Известия Саратовского университета. Сер.: Математика. Механика. Информатика. 2008. Т. 8, Bып. 1. С. 36–42. DOI: http://doi.org/10.18500/1816-9791-2008-8-1-36-42.
  2. Куканов Н.И., Великанов П.Г. Температурный изгиб односвязных и многосвязных изотропных пластин непрямым методом граничных элементов // Актуальные проблемы современной науки: Труды 2-го Международного форума. Естественные науки. Ч. 1–3. Самара: Изд-во СГТУ, 2006. С. 178–183.
  3. Великанов П.Г., Куканов Н.И., Халитова Д.М. Нелинейное деформирование цилиндрической панели ступенчато-переменной жесткости на упругом основании методом граничных элементов // Актуальные проблемы механики сплошной среды. Ижевск: ООО "Принт" 2020. С. 111–115.
  4. Артюхин Ю.П., Грибов А.П. Решение задач нелинейного деформирования пластин и пологих оболочек методом граничных элементов. Казань: Фэн, 2002. 199 с. URL: https://booksee.org/book/438702.
  5. Шевченко В.П. Интегральные преобразования в теории пластин и оболочек. Донецк: Донецкий государственный университет, 1977. 115 с.
  6. Коренев Б.Г. Некоторые задачи теории упругости и теплопроводности, решаемые в бесселевых функциях. Москва: Физматгиз. 1960. 458 с.
  7. Панич О.И. О потенциалах полигармонического уравнения четвертого порядка // Мат. сборник. 1960. Вып. 3. Т. 50. С. 335–354. URL: http://mi.mathnet.ru/msb4795.

Дополнительные файлы

Доп. файлы
Действие
1. JATS XML

© Великанов П.Г., Куканов Н.И., Халитова Д.М., 2021

Creative Commons License
Эта статья доступна по лицензии Creative Commons Attribution 4.0 International License.

Данный сайт использует cookie-файлы

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

О куки-файлах