%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This file is part of NS2DDV.                                                      %
%                                                                                   %
% Copyright(C) 2011-2018    C. Calgaro  (caterina.calgaro@math.univ-lille1.fr)      %
%                           E. Creusé   (emmanuel.creuse@math.univ-lille1.fr)       %
%                           T. Goudon   (thierry.goudon@inria.fr)                   %
%                           A. Mouton   (alexandre.mouton@math.univ-lille1.fr)      %
%                                                                                   %
% NS2DDV is free software: you can redistribute it and/or modify it under the terms %
% of the GNU General Public License as published by the Free Software Foundation,   %
% either version 3 of the License, or (at your option) any later version.           %
%                                                                                   %
% NS2DDV is distributed in the hope that it will be useful, but WITHOUT ANY         %
% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A   %
% PARTICULAR PURPOSE. See the GNU General Public License for more details.          %
%                                                                                   %
% You should have received a copy of the GNU General Public License along with      %
% NS2DDV. If not, see <http://www.gnu.org/licenses/>.                               %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function[p,t] = structmesh_circle_2d(PARAMETERS)

	% INPUT
	% PARAMETERS	List of simulation parameters
	%
	% OUTPUT
	% p 			P1 nodes
	% t				Triangles described by P1 nodes

	nt = 6*PARAMETERS.MESH.NBSEG_C^2;
	np = 3*PARAMETERS.MESH.NBSEG_C*(PARAMETERS.MESH.NBSEG_C+1);     % number of nodes

	t = zeros(3,nt);
	p = zeros(2,np);

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	%%%%                  mesh connectivities array building               %%%%
	ie = 0;
	js = 0;

	for iseg = 1:PARAMETERS.MESH.NBSEG_C
    
		js = js+1;
		jsold = js;
		if iseg==1
			init_is = 1;
		else
			init_is = init_is + 6*(iseg-1);
		end
		is = init_is;
    
		for j=1:6

			for i = 1:iseg-1
            
				is = is+1;
				ie = ie+1;

				t(1,ie) = is;
				t(2,ie) = is+1;
				t(3,ie) = js;

				is = is+1;
				ie = ie+1;

				t(1,ie) = js;
				t(2,ie) = is;
				t(3,ie) = js+1;
            
				js = js+1;
				is = is-1;
            
			end

			is = is+1;
			ie = ie+1;
        
			t(1,ie) = is;
			t(2,ie) = is+1;
			t(3,ie) = js;

		end
    
		if iseg==1
			t(2,ie) = jsold+1; 
		else
			t(3,ie-1)= jsold;
			t(3,ie  )= jsold;
			t(2,ie  )= js;
			js = js-1;
		end
	end

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	%%%%                  Nodes coordonates array building                 %%%%

	is = 1;
	p(1,is) = PARAMETERS.DOMAIN.XCENTER;
	p(2,is) = PARAMETERS.DOMAIN.YCENTER;

	for iseg = 1:PARAMETERS.MESH.NBSEG_C
    
		radius = DOMAIN.RADIUS*iseg/PARAMETERS.MESH.NBSEG_C;
		theta = 2.*pi/(6.*iseg);
    
		for i = 1:6*iseg
			is      = is + 1;
			p(1,is) = radius * cos(theta*(i-1));
			p(2,is) = radius * sin(theta*(i-1));
		end
    
	end

end