Программирование МКЭ в МATLAB. Даутов P.З. - 23 стр.

§ 2. Определение геометрии области, построение P1 сеток                                                         23


i f ∼isempty ( s ) ,
   f o r i s =1:11 % l i n e segments
       i=f i n d ( bs==i s ) ;
       i f ∼isempty ( i )
          x ( i )=xy ( 1 , i s )+( xy ( 2 , i s )−xy ( 1 , i s ) ) * s ( i ) ;
          y ( i )=xy ( 3 , i s )+( xy ( 4 , i s )−xy ( 3 , i s ) ) * s ( i ) ;
       end
   end
   f o r i s =12:15 % c i r c l e segments
       i=f i n d ( bs==i s ) ;
       i f ∼isempty ( i )
          xc =1; yc=−0. 5 ; r c=0 . 2 5 ;
          x ( i )= r c * c o s ( s ( i ))+ xc ;
          y ( i )= r c * s i n ( s ( i ))+ yc ;
       end
   end
end

В начальных комментариях указаны способы вызова функции (эти
способы вызова используются функциями pde toolbox; они реализо-
ваны в выделенном стандартном коде, который не нужно менять).
Определение матрицы d является обязательным: в нем должно быть
столько столбцов, сколько сегментов образует геометрию области. В
нашем случае 15. Для каждого сегмента столбец содержит последова-
тельно значение параметров s0 , s1 , метку областей слева и справа от
сегмента. В матрице xy определены данные для прямолинейных сег-
ментов: 1, 2 (3, 4) строки определяют x (y) координаты начальной и
конечной точки сегмента соответственно. Далее, после стандартного
кода, определяется каждый сегмент в параметрической форме.
    Аналогично можно определить достаточно сложные по форме об-
ласти. Тем не менее, необходимо сделать одно замечание по поводу
параметризации сегментов, поскольку одну и ту же кривую можно
параметризовать множеством способов. В отношении построения се-
ток наиболее удачным является выбор в качестве s параметра длины
дуги сегмента (или пропорционального ему параметра, как например,
выбор полярного угла в случае окружности). Однако такой выбор не
всегда осуществим практически. В этом случае можно поступить сле-
дующим образом.
    Рассмотрим область в виде эллипса с центром в начале координат
и полуосями a = 1, b = 0.1. Напишем файл геометрии gell0 как и выше
и построим сетку. Получим рис. 4:
[ p , e , t ]= i n i t m e s h ( ' g e l l 0 ' , 'Hmax' , 0 . 1 ) ; pdemesh ( p , e , t ) , a x i s e q u a l