Предлагается разновидность коллокационного метода граничных элементов с кубической скоростью сходимости, позволяющего получить решения начально-краевых задач с граничными условиями первого, второго и третьего рода для уравнения dtu = a2A2u - pu с постоянными a,p > 0 в плоской пространственной области при нулевом начальном условии. Для того чтобы иметь возможность доказать сходимость метода с указанной скоростью, аппроксимация интегралов на сингулярных и околосингулярных граничных элементах осуществляется на основе аналитического интегрирования по расстоянию между точками границы. Такая аппроксимация практически и теоретически осуществима для любой аналитически заданной границы класса C5.
On solving plane problems of non-stationary heat conduction by the collocation boundary element method.pdf В настоящей работе предлагается полностью обоснованный коллокационный метод граничных элементов (КМГЭ) [1, с. 21], позволяющий получить численные решения внутренних и внешних начально-краевых задач (НКЗ) с граничными условиями первого, второго и третьего рода для уравнения dtu = а2Д2u - pu с постоянными a, p > 0 в плоской пространственной области при нулевом начальном условии. Решения ищутся в виде потенциала двойного слоя для задачи Дирихле и в виде потенциала простого слоя для задач Неймана и Робина с неизвестными функциями плотности, определяемыми из граничных интегральных уравнений (ГИУ) второго рода. Ранее обоснование коллокационных схем для решения ГИУ второго рода, возникающих в задачах нестационарной теплопроводности, выполнялось в работах [2-4]. В этих работах исследовалась задача Неймана и решение находилось в виде потенциала простого слоя. В работах [2, 4] доказательство сходимости метода было сделано на границах класса гладкости Cш, а в работе [3] -на негладких поверхностях, удовлетворяющих условию типа Липшица. В работе [4] разработана методика численного решения ГИУ с сингулярной правой частью на основе соответствующей замены переменной. Обоснование КМГЭ для плоской задачи Дирихле дано в работе [5], при этом используется ГИУ первого рода и решение ищется в виде потенциала простого слоя. Во всех перечисленных работах, а также других работах, посвященных обоснованию КМГЭ, вопрос аппроксимации интегралов на пространственновременных граничных элементах (ГЭ), которые образуются в результате кусочнополиномиальной интерполяции функции плотности, считается чисто вычислительным и выносится за рамки доказательства сходимости методов, которое предполагает точное вычисление таких интегралов. При этом, как правило, в первую очередь осуществляется аналитическое интегрирование по временной переменной т (при p = 0 это возможно). Затем выполняется интегрирование по длине дуги s (в плоском случае), причем в случае интегрального оператора двойного слоя интегралы вычисляются с помощью простых квадратурных формул Гаусса (ПКФГ) [1, с. 79; 4]. В случае интегрального оператора простого слоя на сингулярных ГЭ выделяется логарифмическая особенность, которая для простых геометрий интегрируется аналитически, а все остальные интегралы по s также вычисляются с помощью ПКФГ [5; 1, с. 172] (сингулярными называются ГЭ, содержащие точку с парой значений т = 0, r = 0; r - расстояние от граничной точки, в которой вычисляется интеграл как функция от параметра, до текущей граничной точки интегрирования). При таком подходе не учитывается неограниченное возрастание производных подынтегральной функции на определенной части околосингулярных ГЭ при измельчении сетки, а также неограниченность таких производных на сингулярных ГЭ в случае интегрального оператора двойного слоя, что не позволяет в принципе применять ПКФГ без снижения порядка аппроксимации схемы в целом. Вместе с тем применение ПКФГ дает хорошие результаты, если точки, в которых вычисляются интегралы на ГЭ как функции от параметра, расположены на границе между ГЭ. Не имея возможности ни обосновать применение ПКФГ теоретически, ни опровергнуть практически, автор предлагает для вычисления интегралов по длине дуги s на сингулярных и околосингулярных ГЭ аппроксимацию третьего порядка относительно шага по длине дуги, использующую осуществимое для любого аналитически заданного контура точное интегрирование по расстоянию r между точками границы. При таком интегрировании роль весовых функций играют функции переменной r, порожденные фундаментальным решением уравнения теплопроводности, а остальная часть подынтегральной функции аппроксимируется с помощью квадратичной интерполяции по г. Так вычисляются интегралы на сингулярных ГЭ, а также на околосингулярных ГЭ в некоторой фиксированной по длине дуги области, прилегающей к сингулярному ГЭ, а на остальных ГЭ интегралы по s вычисляются с помощью ПКФГ. Первоначальное интегрирование по т проводится аналогично: множитель e~p т аппроксимируется с помощью кусочноквадратичной интерполяции, и тогда интегралы вычисляются точно. Стоит отметить, что здесь при решении ГИУ не осуществляется кусочно-полиномиальная интерполяция функции плотности по времени, как это делается обычно в КМГЭ, а выполняется кусочно-квадратичная интерполяция временной С0-полугруппы, через которую выражаются ядра интегральных операторов. В целом получена схема с кубической скоростью сходимости относительно шагов по времени и длине дуги. Доказана сходимость в равномерной операторной топологии как сеточных операторов, аппроксимирующих разрешающие операторы ГИУ, так и сеточных функционалов, определяющих приближенное решение НКЗ в любой точке пространственно-временной области. Доказана устойчивость приближенных решений ГИУ и НКЗ к возмущениям граничных функций. Полученные результаты справедливы для границы с гладкостью C5. Приведены результаты вычислительных экспериментов, подтверждающие кубическую сходимость приближенных решений всех трех НКЗ в круговой пространственной области. В заключение отметим, что в данном КМГЭ используется равномерная временная сетка. Благодаря этому разрешающие ГИУ сеточные операторы экономно вычисляются, так как имеют вид полиномов, образованных степенями оператора правого сдвига на временной сетке. Матричные коэффициенты таких полиномов реккурентно находятся с помощью матричных коэффициентов прямого сеточного оператора ГИУ, имеющего аналогичный вид. Предварительные замечания Пусть O+ - плоская открытая ограниченная односвязная область, и O- = R2\ O+ (R = (-да, +да)). Кроме того, пусть dO, граница области O+ , является кривой класса гладкости C2, если не оговорено особо. Рассмотрим четыре краевые задачи (внутренние и внешние при i = 1,2): а2 Д2u± - pu± = Bu± (x = (x1, x2) e O± ), u± = w± (x e dO ), dnu± - nu± = w± (x e dO ), (1) где u± (x) и w± (x) - векторные функции со значениями в пространстве L2 = L2 (It ) (IT = [0, T]), заданные на множествах Q± и 5Q соответственно (все вводимые здесь пространства функций считаем комплексными); n - нормаль к кривой 30. в точке x e 30, направленная внутрь области Q+; p > 0, Д2 =d2XiX 1 +дХ2 (непрерывность и дифференцируемость векторных функций предполагается здесь в норме пространства их значений, в данном случае - L2 ), a > 0 - коэффициент температуропроводности, ц> 0 - коэффициент теплообмена на боковой поверхности цилиндра (p, a, n - постоянные); B - замкнутый оператор в L2: (Bf) (t) = f '(t), заданный на множестве D(B) абсолютно непрерывных на промежутке IT функций f (t) e L2, таких, что f (0) = 0 . Заметим, что в настоящей работе рассматриваются только линейные пространства и операторы. Пусть C(O') и Ck (O') - пространства непрерывных и к раз непрерывно дифференцируемых на некотором множестве O' с R2 векторных функций со значениями в L2. Авторами [6, 7] доказана однозначная разрешимость задач (1) в классе C(0±) n C2(0±) при любых w± e C(dO). Решения имеют вид векторных потенциалов - криволинейных интегралов первого рода: u± (x) = Gj (x) v± , u± (x) = G0 (x) v± (x e0± ), (2a) где функции v± e C(dO) находятся из соответствующих ГИУ: (G± v±) (x) = w± (x) (x edO , i = 1,2); G± = ±2-1 + Gj, G± = +2-1 + G2-n G0; (2b) Gi (x) f = (G,. f)(x) = f Kt (x, x ') f (x ') ds' (f e C(dO), i = 0,2), dO Ki (x, x') (x ф x') - ограниченные операторы в пространстве L2, определяемые равенствами: K(x,x')f gi (x,x',t)e P%U(x)fdx (f e L2 ), g0(x,x',t) = ao(r,t) , IT gi (x, x', t) = ai (r, t) bi (x, x') (i = 1,2), a0(r, t) = a(r, t) , al(r,t) = a2(r,t) = -r dra(r,t) , bt(x,x') = dn lnr-1. Здесь a(r,t) = (4лт) 'exp|^-r2J(4a2T)J, r = |x-x'|; nl и n2 = n - нормали к кривой дО, проходящие через точки x ' и x соответственно и направленные внутрь области О+; дифференцирования д и д осуществляются соответственно по точкам x ' и x. Операторы U(т) образуют С0 -полугруппу правых сдвигов, порождаемую оператором B: (U(t) f )(t) = f (t -t) при т< t, (U (t) f )(t) = 0 при T>t, Bf = lim t-'( f - U (t) f) (f e D( B)). Заметим, что т^+0 ||U (t)|| = 1 при т< T , U(t) = O при т> T (O - нулевой оператор). Введем в рассмотрение параметрические уравнения кривой дО: xj = Xj(s), x2 = X2 (s). Параметр s по модулю равен длине дуги, откладываемой от некоторой фиксированной точки и заканчивающейся в точке x(s) = ( X1(s), X2(s)), причем s > 0, если дуга откладывается по часовой стрелке, и s < 0, если против. Функции Xj(s), X2(s), периодические с периодом 2S (S- половина длины дО), осуществляют взаимнооднозначное отображение множества IS =[-S, S) на множество дО. Условимся далее писать дО e Ск, если функции Xi (s) (i = 1,2) имеют непрерывные производные на замкнутом множестве IS до порядка к включительно, причем X(l) (-S + 0) = X(l) (S - 0) (l = 1, к). Введем в пространстве С (дО) норму: ||f| С (дО) = max|| f (x)L . Обозначим ( ) xeдQ 2 через Ск (дО) (к e N = {1,2,...}, С0(дО) = С (дО)) банаховы пространства, состоящие из функций f e С (дО), имеющих непрерывные на множестве дО производные f(l): f(l\s) = dl f (x(s) )/dsl (s e IS , l = 1, к), с нормой (l)ll 1С (дО) 1П0гХ 1С (дО) Обозначим через Hn (n e N) гильбертово пространство функций f e L2: -i1/2 Bm f e L2 (m = 1,n), с нормой (H0 = L2). Введем в b и llHn m=01| 2 _ рассмотрение банаховы пространства Ckn (дО) (к e Z+ = {0,1,...} , n e N, С0 (дО) = Ск (дО)), состоящие из элементов f e Ск (дО), таких, что f (x) e Hn при x e дО и Bm f e Ск (дО) (m = 1, n), с нормой С(дО)" l-h -U" . Будем рассматривать также банаховы пространства Скпп1 (дО)=С^ (дО) пС„+т (дО) (к,neZ+ , meN, С^(дО)=Скп (дО)). с нормой С (дО) +l1Л\С„+т (дО) \\Ск„т (дО) = Условимся оператор A, отображающий банахово пространство B в банахово пространство C, обозначать как A [ B ^ C ], а если C = B, то A [ B ]. В силу теорем 2 и 3 [8] имеет место утверждение: Теорема 1. Пусть dQ е Ck + 2. Тогда операторы G± [ Cknm (dQ) ] всюду определены, ограничены и ограниченно обратимы (k,m, n е Z+). Введем в рассмотрение банахово пространство C(Е) непрерывных на множестве Е = dQx IT скалярных функций f (x, t) с нормой ||f|C(Е) = maxjf (х, t)|. (е) (x,t)eE Введем также в рассмотрение банаховы пространства Ck(Е) и Cn(Е) (k е Z+, n е Z+, C0(Е) = C0(E) = C(Е)) функций f е C(E), имеющих непрерывные на множестве IS xIT скалярные производные (dlsf)(s,t) = dlsf (x(s),t) (l = 0,k) и dm f (X,t) (m = 0,n) соответственно с нормами i-i. = max d\ и llCk (Е) C (Е) l=0,k . Наконец, будем рассматривать банаховы пространства m = max d t "Cn (Е) C (Е) C^ (Е) = Ck (Е) n Cn (Е) с нормой (k, n е Z+). С уче- C (Е) +11 C (Е) "Cn (Е) том соотношений lV2 |f (x, t)| = j(Bf )(x, t’)dt’< t\\(Bf )(x, t')\2 dt’ 0; p = -r , если ст< 0 . Введем в рассмотрение функции уг- (s, s') (i = 0,3), заданные на множестве IS x IS при s ' Ф s равенствами уг- = фг/ст2 (i = 0,2) и ф3 =ф3/ст, где ф0(s, s')= r 2 = [ (s ') - Xi (s)]2 + [X2 (s ') - X2 (s)]2 , Ф1 (s, s ') = r dni r = x2 (s' ) [Xi (s ') - Xi (s)] - x'i (s ') [X2 (s ') - X2 (s)] , Ф2 (s, s ') = r dni r = x2 (s) [Xi (s) - Xi (s ')] - x’i (s) [X2 (s) - X2 (s ')] , Ф3(s, s ') = r dsr = xi(s ') [xi(s' ) - xi(s)]+ x2(s' ) [x2(s ') - x2(s)], а при s = s равенствами Ф0 (s, s) = ф3 (s, s) = i , (s, s) = ф2 (s, s) = 2-i [x[(s) x2 (s) - x2 (s) x| (s)]. Кроме того, введем в рассмотрение функции 5(s,s') = (ds'р)-i ^/ф^/ф3 , Ь,(s,s ') = bt(x(s),x(s')) = -\Vt/\V0 (i = i,2). Лемма [8]. Пусть I - замкнутый интервал на вещественной оси. Предположим, что некоторая вещественная функция f (z, Q имеет на множестве I х I непрерывные производные д\д^ f (i = 0, m , j = 0,m’), причем m < m’ и djf \^=z = 0 при z e I, j = 0, q -1, где q = ml - m. Тогда функция h(z, Q , заданная при Z ^ z равенством h(z,Z) = f/(Z-z)q , а при Z = z равенством h(z,z) = dq^f \z=z/q!, имеет на множестве I х I непрерывные производные д\dljh при i = 0, m - j , j = 0, m . Теорема 2. Пусть dQe Cn+2 (n e Z+). Тогда существуют непрерывные на множестве IS хIS производные д]s,bi (j = 0,n , i = 1,2 ). Кроме того, для любого M > 1 существует Е > 0, такое, что при (s, ст) e IS х ЬЕ (ЬЕ = [-Е, Е]) функция 5 ограничена: 1 1 существует Е > 0 , такое, что 1 JcrfM при (s, ст) e Is х IE и получаем утверждения теоремы на основе представлений дs'5 = Fj /(фз +V0-V2 ) ( ст) e IS х IE X дsb = Gi, j I+1 ((s, s') e Is х Is , i = 12X где функции Fj и Gi j суть линейные комбинации произведений некоторых степеней производных дk ф3, д1 -у0 (k, l = 0, j) и дk ф, д1 -у0 (k, l = 0, j ) соответственно. Теорема доказана. Следствие 1. Пусть дО e Cn+2 (n e Z+). Тогда функция ps (ст) = p(s, s + ст) при любых фиксированных s e Is и M > 1 диффеоморфно с гладкостью Cn+2 отображает множество IE на множество ps (IE) = [ps (-Е), ps (Е)]. Функции b0(s, р) = 5(s, s + ст s (р)), b (s, p) = 5(s, s + ст s (p)) (s, s + .^ (p)) (s e Is , i = 1,2), где us (p) - функция, обратная к функции ps (ст), имеют непрерывные на множестве Is х ps (ЬЕ) производные д1pbi (j = 0, n , i = 0,2). Обозначим через Лт (z) и Лm (z) (z е [a, b], m = 0,2) квадратичные функции - интерполяционные многочлены Лагранжа, определенные на промежутке [a,b]: 2 z _ z ■ Лm (z) ■ П --- > zm ■ z + qmhz '; }=0(j*m) zm Zj Л m (z) = П - _ z , Z m = z + qmhz ( m = 0,2 )j=0(j*m) Zm Zi Здесь hz ■ 2~l(b _a), z ■ 2_ (a + b); q0 ■ _1, q1 ■ 0, q2 ■ 1; q0 ■ _-n/I/2 , q1 ■ 0, q2 ■ VI/2 [9, с. 92]. Пусть f (z) - трижды непрерывно дифференцируемая на промежутке [a, b] функция со значениями в произвольном банаховом пространстве B . Тогда для функций f1(z) = Z m=0 f (zm (z) и f2 (z) = S m=0 f (zm )Лm (z) ( z е [a, b] ) имеют место оценки: lf _ -4[ab] < C.V ?f(1)(z)dz (z0 _ z2 )(z0 _ z1 )• 4 + 2z _(z" +z1) | f(1>(grfC, (z2 _ z0 )) z2 _ z1 ) -) f (2)(C))z0 _Z)dZ z1 К 1 z0 ))z2 1 )z0 _ z2 ))z0 f/^ z) - f f (2)(Q)z2_Z)d. zr* Аналогичные формулы имеют место и для функции (z). Вследствие этих фор мул имеют место оценки ||/П) 0; s't = max , -E}, s” = min ,-E}, если sl < 0, а чис ло E > 0 выбрано в соответствии с теоремой 2. В интегралах J'imkl(т) на основании следствия 1 сделаем замену переменной ст = СТ (р): р (si+i) ^ _ _ Jim,k,l (Т) = I ai (P, Т) bi,m,k (Р) dР( h m = 0,2, k, l = -L, L - II Pk(s!) bi,m,k (P) = b (sk , Р)Лm (k + CTk (P) ) , Pk (ct) = p(sk, sk + ст), стк (p) - функция, обратная к функции pk (а). Введем в рассмотрение соответствующие интегралы, аппроксимирующие Jimkl(т): Rt(s'i+i) л _ _ Ji,m,k ,l (т) = { ai (P, т) bi,m,k (P) d P( i, m = 0,2,k, l =-4 L - 1 ), Pk(s/) 2 _ b,m,k (P) = Z b>,m,k (Rk,i + qmhk,l )Лm' (P) ( P6 [Pk (sl), Pk (sI+l)] ), m '=0 hk,i = 2-1 [Pk (si'+i) - Pk (si)], Pk,i = 2-1 [Pk (si) + Pk (sM)]. Интегралы J"mki(т) аппроксимируем с помощью ПКФГ с у узлами: Ji,m,k,l (т) = 2-1 hs Z Wj g,m,k (Sl + 2-1 hs Zj , т) , Sl = 2-1 (sl + Sl+1 ) ( h m = 0, 2 , j=1 k,l = -L, L -1), при этом Zj - корни многочлена P,f (z) = [y !/(2y)!](d(dzY)z2 -1) на промежутке [- 1; 1] [9, с. 258]; для весовых коэффициентов Wj выполняется равенство Wj = 2 (Wj > 0) [9, с. 255]. Операторы Gi, Gi n (n = 0, N -1), в которых интегралы Jimk{ (т) заменены выражениями J’i m k,l (т) и J”m k,l (т), обозначим через Gi, G'in и G’, G” n соответственно. В силу следствия 1, теоремы 2 и неравенства r > cr Е , имеющего место при |ст|>Е, при указанной гладкости границы могут быть определены константы cij = max \d]pbi(s,p)| (5Qe Cn+2), С",} = _sup l^g*(s,s + ст,т)| (s,p)eIs xps (Ie) s6ls ,ст>Е,т>0 (dQe Cn) (j = , i = 0,2). Учитывая неравенства (4b), (6) и hk l < 2-1 hs, при условиях 5Q 6 C5 и f 6 C2(5Q) при любых i = 0,2, k =-L, L -1, l = -L/ 2,L/2-1 и т>0 имеем оценку Pit (s2l + 2-k ) IK Ji ,m,k .2l+l'- k (т) - Ji ,m,k ,2l+l '-к (т) )f2l + m < =01'=0 L2 ft (s2l + 2-к ) (,max, J5P[*i (sk . P) f (sk + CTk (P) )]|| L f ai (P.т) d P< PePk ([ - k .*21+ 2-k ])) IIL2 J < 8-1 h cffl Pk (s2l-k ) c h* ii / iu2 s II./ IIC2(3Q) f ai (P. t) dP . Pk (s2l-k ) ci - 8 Cffl [ci .3 СЛ + (3ci.2 c0,0 + 3ci .1 c0.1 + Ci .0 c0.2 )СЛ + 3 (Ci .1 c0.0 + Ci .0 c0.1 c0.0 ) СЛ ] . (12a) где f = f (x{). f - PL PL f. Если 5Q e C2 . то при любых i = 0.2. k. l = -L. L -1 и т> 0 имеем оценку 2 1 < C"hs2y+1H f| C 2(SO)> IK Ji.m.k .2l+Г - k (т) - Ji .m.k .2l+l '-k (т)) m=01'=0 -^2 (12b) , „ - (Y !)4 [
Бреббия К., Теллес Ж., Вроубел Л. Методы граничных элементов. М.: Мир, 1987. 524 с.
Onishi K. Convergence in the boundary element method for heat equation // Teaching for Robust Understanding of Mathematics. 1981. V. 17. P. 213-225.
Costabel M., Onishi K., Wend land W.L. A boundary element collocation method for the Neumann problem of the heat equation // Inverse and Ill-Posed Problems (H.W. Engl and C.W. Groetsch, ed.). Boston: Academic Press, 1987. P. 369-384.
Hongtao Y. On the convergence of boundary element methods for initial-Neumann problems for the heat equation // Mathematics of Computation. 1999. V. 68. No. 226. P. 547-557.
Hamina M., Saranen J. On the spline collocation method for the single layer heat operator equation // Mathematics of Computation. 1994. V. 62. No. 205. P. 41-64.
Иванов Д.Ю. Решение двумерных краевых задач, соответствующих начально-краевым задачам диффузии на прямом цилиндре // Дифференц. уравнения. 2010. Т. 46. № 8. С. 1094-1103.
Иванов Д.Ю., Дзержинский Р.И. Решение задач Робена для двумерных дифференциально-операторных уравнений, описывающих теплопроводность в прямом цилиндре // Научно-технический вестник Поволжья. 2016. № 1. С. 15-17.
Иванов Д.Ю. Устойчивая разрешимость в пространствах дифференцируемых функций некоторых двумерных интегральных уравнений теплопроводности с операторно-полугрупповым ядром // Вестник Томского государственного университета. Математика и механика. 2015. № 6 (38). С. 33-45.
Березин И.С., Жидков Н.П. Методы вычислений. Т. 1. М.: ГИФМЛ, 1962. 464 с.
Фихтенгольц Г.М. Курс дифференциального и интегрального исчисления. Т. 2. М.: ФИЗМАТЛИТ, 2001. 810 с.
Иванов Д.Ю. Вычисление операторов, разрешающих задачи теплопроводности в прямых цилиндрах, с использованием полугрупповой симметрии // Известия Московского государственного технического университета МАМИ. 2014. Т. 4. № 4(22). С. 26-38.