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

UptoLike

Составители: 

§ 4. Сопряженная кодировка сетки 31
[ ep , , j ]=unique ( s o r t ( ep , 2 ) , ' rows' ) ; % unique mesh edge s
ep=ep ' ; j=j ' ; nep=s i z e ( ep , 2 ) ; % ep i s s o r t ed : 1 row and each column
% s e t te
mt = 1 : nt ;
te = [ j (mt ) ; j (mt+nt ) ; j (mt+2
*
nt ) ] ; % el eme nt s ed ges
te = [ t e ; t ( 4 , : ) ] ;
% s e t et = 6 ,7 rows o f ee ( f i r s t attempt )
et =z e r o s (2 , nep ) ; % ed ge2element c o n n e c t i v i t y matrix
f o r k=1: nt
f o r j =1:3
i t=te ( j , k ) ;
i f et ( 1 , i t )==0, et ( 1 , i t )=k ;
e l s e et ( 2 , i t )=k ; end
end
end
% s o r t f i r s t 2 rows o f e as i n ep
f o r j =1: ne
i f e ( 1 , j )>e ( 2 , j ) ;
c=e ( 1 , j ) ; e ( 1 , j )=e ( 2 , j ) ; e ( 2 , j )=c ;
c=e ( 3 , j ) ; e ( 3 , j )=e ( 4 , j ) ; e ( 4 , j )=c ;
c=e ( 6 , j ) ; e ( 6 , j )=e ( 7 , j ) ; e ( 7 , j )=c ;
end
end
[ , ib , i t ] = i n t e r s e c t ( ep ' , e ( 1 : 2 , : ) ' , ' rows' ) ;
% s e t 3 : 5 rows o f ee
ee =[ep ; zer o s ( 5 , nep ) ] ;
ee ( 3 : 5 , ib )=e ( 3 : 5 , i t ) ;
% s e t 6 ,7 rows o f ee
f o r i i =1: nep % i i =edge
t1=e t (1 , i i ) ; t2=e t ( 2 , i i ) ; % n ei gh bo ur s t r i a n g l e s
i=ep ( 1 , i i ) ; j=ep ( 2 , i i ) ; k=s e t d i f f ( t t ( t1 , : ) , [ i j ] ) ; % t1 v or te x i n d i c e s
% o r i e n t e d are a o f t r i a n g l e ( i , j , k )
d12 = p ( : , j )p ( : , i ) ;
d13 = p ( : , k)p ( : , i ) ;
A = d12 ( 1 , : ) .
*
d13 (2 ,: ) d12 ( 2 , : ) .
*
d13 ( 1 , : ) ; % A/2=area
i f A>0, ee (6 , i i )= t1 ; ee ( 7 , i i )= t2 ;
e l s e ee ( 6 , i i )=t2 ; ee ( 7 , i i )= t1 ; end
end
Для примера выполним следующие команды:
[ p , e , t ]= i nit me sh (' l sha peg ' , 'hmax' , i n f ) ;
[ ee , te ]=adjmeshP1 (p , e , t ) ;
plotadjmeshP1 (p , e , t , ee , 1 2 ) ;
где функция рисования имеет следующий вид:
§ 4. Сопряженная кодировка сетки                                                                                     31


[ ep ,∼ , j ]= unique ( s o r t ( ep , 2 ) , 'rows' ) ;       % unique mesh e d g e s
ep=ep ' ; j=j ' ; nep=s i z e ( ep , 2 ) ; % ep i s s o r t e d : 1 row and each column

% set te
mt = 1 : nt ;
t e = [ j (mt ) ; j (mt+nt ) ; j (mt+2* nt ) ] ; % e l e m e n t s e d g e s
te = [ te ; t ( 4 , : ) ] ;

% s e t e t = 6 , 7 rows o f e e ( f i r s t attempt )
e t =z e r o s ( 2 , nep ) ; % e d g e 2 e l e m e n t c o n n e c t i v i t y matrix
f o r k=1: nt
    f o r j =1:3
        i t=t e ( j , k ) ;
        if      e t ( 1 , i t )==0, e t ( 1 , i t )=k ;
        e l s e e t ( 2 , i t )=k ; end
    end
end

% s o r t f i r s t 2 rows o f e a s i n ep
f o r j =1: ne
    i f e ( 1 , j )>e ( 2 , j ) ;
       c=e ( 1 , j ) ; e ( 1 , j )=e ( 2 , j ) ; e ( 2 , j )=c ;
       c=e ( 3 , j ) ; e ( 3 , j )=e ( 4 , j ) ; e ( 4 , j )=c ;
       c=e ( 6 , j ) ; e ( 6 , j )=e ( 7 , j ) ; e ( 7 , j )=c ;
    end
end
[∼ , ib , i t ] = i n t e r s e c t ( ep ' , e ( 1 : 2 , : ) ' , 'rows' ) ;

% s e t 3 : 5 rows o f e e
e e =[ ep ; z e r o s ( 5 , nep ) ] ;
e e ( 3 : 5 , i b )=e ( 3 : 5 , i t ) ;

% s e t 6 , 7 rows o f e e
f o r i i =1: nep % i i =edge
    t 1=e t ( 1 , i i ) ; t 2=e t ( 2 , i i ) ; % n e i g h b o u r s t r i a n g l e s
    i=ep ( 1 , i i ) ; j=ep ( 2 , i i ) ; k= s e t d i f f ( t t ( t1 , : ) , [ i j ] ) ; % t 1 v o r t e x i n d i c e s
   % oriented area of t r i a n g l e ( i , j , k)
    d12 = p ( : , j )−p ( : , i ) ;
    d13 = p ( : , k)−p ( : , i ) ;
   A = d12 ( 1 , : ) . * d13 ( 2 , : ) − d12 ( 2 , : ) . * d13 ( 1 , : ) ; % A/2= a r e a
    i f A>0, e e ( 6 , i i )= t 1 ; e e ( 7 , i i )= t 2 ;
    else          e e ( 6 , i i )= t 2 ; e e ( 7 , i i )= t 1 ; end
end

Для примера выполним следующие команды:
[ p , e , t ]= i n i t m e s h ( ' l s h a p e g ' , 'hmax' , i n f ) ;
[ ee , t e ]= adjmeshP1 ( p , e , t ) ;
plotadjmeshP1 ( p , e , t , ee , 1 2 ) ;

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