19.4.1999

>> pdecirc(0,0,2)

Session saved by means of GUI file-menu This session can be edited and run again.

After adjusting the boundary conditions we have:

>> [p,e,t]=initmesh(g,'Hmax',5,'jiggle','off');
>> pdeplotn(p,e,t)

Look at the boundary edges first. 
e([1,2,5],:)

ans =

     1     3     2     7     5     8     1     6
     2     4     7     5     8     4     6     3
     1     2     3     3     3     3     4     4



Comparing the figures above and below to the values we see (both
geometrically and "algebraically") that
edge 1->2 has boundary label 1
edge 3->4 has boundary label 2
edge 2->7 has boundary label 3
etc.

                  

                  



 Save, edit, run sess2.m and export 
sess2.m again.
>> sess2         
%  see pdesetbd(BOUNDS,TYPE,MODE,COND,COND2,COND3,COND4)
%      PDESETEQ(TYPE,C,A,F,D,TLIST,U0,UT0,RANGE)
>> whos
  Name      Size         Bytes  Class

  a         1x3              6  char array
  b        12x4            384  double array
  c         1x3              6  char array
  d         1x3              6  char array
  f         1x3              6  char array
  g        10x4            320  double array

Grand total is 100 elements using 728 bytes

 Remember the p,e and t 
>> [p,e,t]=initmesh(g,'Hmax',5,'jiggle','off');
>> pdeplotn(p,e,t)
>> p10=p(:,10)

p10 =

    1.2716
    0.8497

>> plot(p10(1),p10(2),'o')

We see that we will get a circle right at number 10, thus indeed the 10th
 point of p is point nr. 10 in the figure (a miracle).

Let's plot the global linear basis function associated to node 10.
We just need to create a vector z s.t. z(10)=1 and z(i)=0 otherwise.
>> z=eye(11,11);
>> z=z(:,10);
>> figure
>> pdesurf(p,t,z)
>> view(-40,80)
>> hold on
>> pdeplotn(p,e,t)


                  

                  

 How to check whether a (set of) node(s) is on the boundary 

The e-matrix has the boundary information. The first two columns 
give the boundary edges (compare to the figure).

>> e12=e(1:2,:)

ans =

     1     3     2     7     5     8     1     6
     2     4     7     5     8     4     6     3


"Ravel" columnwise and transpose:
>> e12=e12(:)'

e12 =
  Columns 1 through 12 
     1     2     3     4     2     7     7     5     5     8     8     4
  Columns 13 through 16 
     1     6     6     3

>> bp=p(:,e12);

>> plot(bp(1,:),bp(2,:),'*')


This is ok, except there are repetitions. Just write (an elegant) function
to remove duplicates (from e12). (There is a standard idiom in APL for this, but Matlab
is not exactly APL.)

 The Matlab function  unique   
As a matter of fact there is a function that does it right away for us:
e12=unique(e12)
bp=p(:,e12)
axis([-.1 2.1 -.1 2.1])  
plot(bp(1,:),bp(2,:),'*')
axis([-.1 2.1 -.1 2.1]) 




There is a function for the boundary testing written by JvP also, but for the 
present purpose the
above is good enough (and gets us used to our friends p, e and t).

Check here! 

Function onboundary.m

Other "course tools"
%function  [isonbnd,index]=onboundary(indp,e)
% 
% indp = indexes of points (=of columns of matrix 'p')
% boudary information matrix 'e'
% 
% Returns the maximum number of points of each column of indp 
% that lie on boundary.
%
% E.g. To test whether some edge of the triangle tr=[p0;p1;p2] 
% lies on boundary we say isonbnd(tr,e)==2. 
% 
% in 'index' we return the column indexes in 'e' matrix for corresponding 
% nodes in 't'; the nodes not on boundary have zero index.
% By J-v-P -97


Here's how the function could be used in our present case.
npts=size(p,2)
[isonbnd,index]=onboundary(1:11,e)
isonbnd==2
p(:,isonbnd==2)
bp=p(:,isonbnd==2)  % Boundary points
plot(bp(1,:),bp(2,:),'*')
axis ([-.1 2.1 -.1 2.1])
figure
pdeplotn(p,e,t)

                  




Continue with linassembly