1 / 6

% SET CVFEM COEFICIENTS----------------- %zero out coefficient for i=1:n %loop on nodes

% SET CVFEM COEFICIENTS----------------- %zero out coefficient for i=1:n %loop on nodes %zero out coefficient ap(i)=0.0; %node point coefficient V(i)=0.0; %volume of CV for j=1:n_s(i)+1 a(i,j)=0.0; %neighbor coefficient end % Set loopcount

haruko
Download Presentation

% SET CVFEM COEFICIENTS----------------- %zero out coefficient for i=1:n %loop on nodes

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. % SET CVFEM COEFICIENTS----------------- %zero out coefficient for i=1:n %loop on nodes %zero out coefficient ap(i)=0.0; %node point coefficient V(i)=0.0; %volume of CV for j=1:n_s(i)+1 a(i,j)=0.0; %neighbor coefficient end % Set loopcount % At internal nodes--recognized by S(i,n_s(i)+1)<>0 % the nodes in the support match the elements in the support % Hence loop over every node in the node i support: loopcount=n_s(i) % At boundary nodes--recognized by S(i,n_s(i)+1)=0 % the nodes in the support are one larger than the elements in the support % Hence only loop over the first n_s-1 nodes: loopcount=n_s(i) loopcount=n_s(i); if S(i,n_s(i)+1)==0 loopcount=n_s(i)-1; end

  2. % Pick out each element in turn from i_th support for j=1:loopcount k1=S(i,j); %The numbers of the nodes in the j_th elemnt are k2=S(i,j+1); % i, k1=S(i,J), and k2=S(i,j+1) xL(1) = x(i); %xL=local x-coordinate of a node within a element xL(2) = x(k1); xL(3) = x(k2); yL(1) = y(i); %yL=local y-coordinate of a node within a element yL(2) = y(k1); yL(3) = y(k2); % Vele= element-volume Vele=(xL(2)*yL(3)-xL(3)*yL(2)-xL(1)*yL(3)+xL(1)*yL(2)... +yL(1)*xL(3)-yL(1)*xL(2))/2; Nx(1)= (yL(2)-yL(3))/(2*Vele); %derivative of shape function wrt x Nx(2)= (yL(3)-yL(1))/(2*Vele); Nx(3)= (yL(1)-yL(2))/(2*Vele); Ny(1)=-(xL(2)-xL(3))/(2*Vele); %derivative of shape function wrt y Ny(2)=-(xL(3)-xL(1))/(2*Vele); Ny(3)=-(xL(1)-xL(2))/(2*Vele);

  3. % Face 1--face area normal components delx= [(xL(1)+xL(2)+xL(3))/3]-[(xL(1)+xL(2))/2]; dely= [(yL(1)+yL(2)+yL(3))/3]-[(yL(1)+yL(2))/2]; % Face 1--Diffusion Coefficients ap(i) = ap(i) + (- Nx(1) * dely + Ny(1) * delx); a(i,j) = a(i,j) + ( Nx(2) * dely - Ny(2) * delx); a(i,j+1) = a(i,j+1) + ( Nx(3) * dely - Ny(3) * delx); % Face 2--face area normal components delx= [(xL(1)+xL(3))/2]-[(xL(1)+xL(2)+xL(3))/3]; dely= [(yL(1)+yL(3))/2]-[(yL(1)+yL(2)+yL(3))/3]; % Face 2--Diffusion Coefficients ap(i) = ap(i) + (- Nx(1) * dely + Ny(1) * delx); a(i,j) = a(i,j) + ( Nx(2) * dely - Ny(2) * delx); a(i,j+1) = a(i,j+1) + ( Nx(3) * dely - Ny(3) * delx); end %end loop on nodes/elements in support % For internal nodes LAST node in region of support is also the FIRST % Hence coefficient for last node needs to be added to first node % and then removed a(i,1)=a(i,1)+a(i,n_s(i)+1); a(i,n_s(i)+1)=0; end %end loop on nodes in domain

  4. The above will generate coefficients in the equation—assuming that there is no flow across segment of Control Vol along boundary i Need This term In a simple treatment we write boundary contribution at node k in form And equation takes the form (1) If boundary is insulated that we set If boundary has fixed value T =Tfix Then we set Where A large value—forcing solution of (1) at node i to return value

  5. %SET BOUNDARY for i=1:n %Default values BC(i)=0; %coefficient part BB(i)=0; %constant part end % Set value = 0 on r = r_outer boundary segment 2 for k=1:n_b(2) BC(B(2,k))=1e18; end %set value =1 on r = r_inner boundary segment 4 for k=1:n_b(4) BC(B(4,k))=1e18; BB(B(4,k))=1e18; end %------------------------------

  6. Test problem steady sate heat conduction in annulus BUT Solved in Cartesian Form in With boundary conditions Using Grid generation program as given in notes –develop a CVFEM Matlab Code (internal coefficients, boundary coefficients, solver) By Tuesday 17 hand complete matlab program and plot that comparing numerical sol With analytical

More Related