
function demo_cubature_polyhedron_01(exampleV,ade,extraction_type)

%--------------------------------------------------------------------------
% Object:
%
% Cubature test over polyhedra (convex or concave, even with odd shapes).
%
% This routine was tested on
%
% Matlab 9.5.0.1298439 (R2018b) Update 7
%
% MacBook Pro (Retina, 13-inch, Early 2015), 2,7 GHz Intel Core i5,
% 16 GB 1867 MHz DDR3
%
% To launch the demo, just digit in Matlab Command window:
%
% >> demo_cubature_polyhedron_02
%
% and "press enter".
% 
% A battery of tests will be performed, on different polyhedra, for the
% methods working with partitioning in tetrahedra or by divergence theorem.
% The key point is the evaluation of monomial moments.
%--------------------------------------------------------------------------
% Warning: 
%
% If ade > 10, then "demo_cubature_polyhedron" can be very time consuming 
% also offering bad compression results.
%--------------------------------------------------------------------------


if nargin < 1
    exampleV=-8:13;
end

if nargin < 2
    ade=6;
end

if nargin < 3
    extraction_type=0;
end

diary on;

for example=exampleV
    test_cubature_polyhedron_01_sub(example,ade,extraction_type);
    fprintf('\n \n \t ................. PAUSE (2 seconds) ................. \n \n');
    pause(2);
end

diary off;







function test_cubature_polyhedron_01_sub(example,ade,extraction_type)

% .............................. settings .................................
if nargin < 1
    example=-1;
end

if nargin < 2
    ade=6; % algebraic degree of exactness
end

tri_method=1; % method to determine the facets: 1: alphaShape 2: convhull

% ............................. main code .................................

% Vertices of the polyhedron.
[vertices,polyhedron_type,alphashape_radius]=provide_vertices(example);

% Facets of the polyhedron.
X=vertices(:,1); Y=vertices(:,2); Z=vertices(:,3);

fprintf('\n \t *** determining facets and triangulation: ');
tic;
switch tri_method
    
    case 1
        % fprintf('\n \t ... USING ALPHA SHAPE ...');
        SHP = alphaShape(X,Y,Z,alphashape_radius);
        facets = boundaryFacets(SHP);
        tri = alphaTriangulation(SHP);
        
    case 2
        % fprintf('\n \t ... USING CONVHULL ...');
        facets = convhull(X,Y,Z);
        tri = delaunay(X,Y,Z);
        
end
toc;


if size(tri,1) == 0
    warning('Not enough points to determine polyhedron');
    return;
end


% Determining rule with total degree "ade"

% Full rule

fprintf('\n \t *** cubature_polyhedron_3Dtriangulation_01 (full):  ');

tic;
[XYZW]=cubature_polyhedron_3Dtriangulation_01(vertices,ade,tri);
toc;

% Full rule + compression

fprintf('\n \t *** cubature_polyhedron_3Dtriangulation_01 (full+comp):  ');

tic;
[XYZW,XYZW_C]=cubature_polyhedron_3Dtriangulation_01(vertices,...
    ade,tri);
toc;


tic;
[XYZW_TC,flag,iteration_stats]=cubature_polyhedron_free(vertices,...
    facets,ade,extraction_type);
fprintf('\n \t *** cubature_polyhedron_free:  ');
toc;

% evaluation of monomial moments up to total degree "ade"

fprintf('\n \t *** test: monomial moments computation: ');
tic;
monomial_moments_XYZW=evaluate_monomial_moments(XYZW,ade);
monomial_moments_XYZW_C=evaluate_monomial_moments(XYZW_C,ade);
monomial_moments_XYZW_TC=evaluate_monomial_moments(XYZW_TC,ade);
toc;


% .............................. STATISTICS ...............................

fprintf('\n');

fprintf('\n \t EXAMPLE: %2.0f ADE: %3.0f',example,ade);
str=strcat('[',polyhedron_type,']'); fprintf('\n \t',str);
dim=(ade+1)*(ade+2)*(ade+3)/6;
fprintf('\n \t VERTICES: %5.0f FACETS: %5.0f TETRAHEDRONS: %6.0f',...
    size(vertices,1), size(facets,1), size(tri,1)); 

err0=norm(monomial_moments_XYZW-monomial_moments_XYZW_C,inf);
err1=norm(monomial_moments_XYZW-monomial_moments_XYZW_TC,inf);

tolmoms=10^(-10);
ifilt=find( (abs(monomial_moments_XYZW) > tolmoms) & ...
    (abs(monomial_moments_XYZW) < 1/tolmoms)) ;

mm_filt=monomial_moments_XYZW(ifilt);
mmC_filt=monomial_moments_XYZW_C(ifilt);
mmTC_filt=monomial_moments_XYZW_TC(ifilt);

Rerr0=norm((mm_filt-mmC_filt)./mm_filt,inf);
Rerr1=norm((mm_filt-mmTC_filt)./mm_filt,inf);

fprintf('\n \t MONOMIAL MOMENTS tri/cmpx (infnorm) AE: %1.1e RE: %1.1e',...
    err0,Rerr0);
fprintf('\n \t MONOMIAL MOMENTS tri/free (infnorm) AE: %1.1e RE: %1.1e',...
    err1,Rerr1);
fprintf('\n');
switch extraction_type
    case 0
        str_et=' Dessole-Marcuzzi-Vianello algorithm';
    case 1
        str_et=' lsqnonneg algorithm';
    case 2
        str_et=' Slawski algorithm';
    case 3
        str_et=' QR algorithm';
end
str_et1=strcat(' Compression by',str_et);
fprintf('\n \t'); disp(str_et1)

for i=1:size(iteration_stats,1)
    fprintf('\n \t EXTRACT: iteration: %2.0f pts: %6.0f in: %6.0f res: %1.1e iters: %3.0f',...
        iteration_stats(i,1),iteration_stats(i,2),iteration_stats(i,3),...
        iteration_stats(i,4),iteration_stats(i,5));
end

fprintf('\n');

% fprintf('\n \t FULL RULE : MONOMIAL MOMENTS ERR. (inf norm): %1.1e',err);
fprintf('\n \t ----------------------------------------------------------');
fprintf('-------- \n \n');


plot(SHP)

function [vertices,polyhedron_type,alphashape_radius]=provide_vertices(example)

alphashape_radius=realmax;

switch example
    
    case -8
        
        polyhedron_type='Concave Polyhedron: many vertices';
        
        vertices=[
            3.477297788921089e-01 8.418478950828894e-01 3.541900202654923e-01
            7.275258801332074e-02 9.752662918392632e-01 3.297929391749036e-01
            5.444395606010217e-01 2.821742380004579e-01 7.270570492677413e-01
            9.653957766977979e-01 9.755712153539947e-01 6.983560482827161e-01
            2.811729756753140e-01 6.021704305592491e-01 2.787669164999301e-01
            8.013582088441695e-02 1.869449155425372e-01 8.671548051229186e-01
            1.458015089554773e-01 4.633294088033622e-01 8.529946605721650e-02
            1.630456215691828e-01 1.241671541812516e-01 6.831448927681127e-02
            5.385779351763819e-01 6.999650618459788e-01 9.682528743769193e-01
            4.689612051040591e-01 7.584145963171395e-01 6.250752166155358e-01
            2.653143384311523e-01 4.797478178312440e-01 7.701803945873731e-01
            6.745387431109657e-01 6.831209378569906e-01 1.453881584285108e-01
            1.245445288320554e-01 7.180495846219421e-01 3.138500378985871e-01
            2.517656984562506e-01 2.127831172843369e-01 6.557461516599750e-04
            5.891697292688555e-01 1.404202343937110e-01 5.594565887876459e-01
            3.510981283107595e-01 4.806924102296263e-01 9.221222383455038e-01
            5.227738592557006e-01 9.994071279610198e-01 4.453218212923998e-01
            3.351707114077863e-01 8.807975766648165e-01 3.492652968108898e-01
            2.353512537843550e-02 3.892564899118944e-02 5.130612542987076e-01
            4.016796381890281e-01 9.486334359463053e-01 8.928140949936025e-02
            5.046834838952802e-01 8.928604398190721e-01 7.637834742008680e-01
            1.295893795789318e-01 2.123250312303010e-01 3.883933120023740e-01
            6.129693985681371e-01 9.819487436278350e-01 1.478423955582254e-01
            9.634558272788004e-01 2.545469919845407e-01 3.850670664619223e-01
            8.965713477692505e-02 6.268564350308548e-01 2.725114083510977e-01
            3.086620967131966e-01 5.187671074096761e-01 2.992473140772639e-01
            5.292138565934588e-01 4.948410460497292e-01 1.153282643043885e-02
            4.474150693200879e-01 2.466037515869339e-01 2.388758920249922e-01
            9.399629070177261e-01 6.375590931464437e-02 4.618937833220929e-01
            7.211043294462249e-01 1.684435235388758e-01 1.408395796809047e-01
            7.378828324064585e-01 6.648762188841058e-01 3.758150551691928e-01
            6.391319850336052e-01 7.060781725166326e-01 4.526165350958380e-01
            3.531449479589577e-02 8.580702129685412e-01 5.268458720123642e-01
            6.416154641968342e-01 2.698257521082746e-01 5.513997158392359e-01
            3.185906420401339e-01 2.246929099482042e-02 3.323731620561032e-01
            1.011017478197296e-01 5.581303275024840e-01 2.149559763769792e-01
            8.098216288263882e-02 1.989752359133103e-01 1.491260723610787e-01
            9.021453571638910e-01 9.202411567085241e-01 5.604570059594651e-01
            2.527356920226418e-01 2.750981001621488e-01 3.286394731326676e-01
            4.949077580412402e-01 2.289603216817687e-01 6.536758968585203e-01
            3.881130254183480e-01 5.181496732431041e-01 7.914050159375723e-01
            2.636263282095280e-01 2.560600808444836e-01 8.108992463142731e-01
            2.631162052449976e-01 1.614860599523993e-02 6.829282685847073e-01
            6.995274292064517e-01 9.764857535311057e-01 2.825646554758101e-01
            5.673144686645981e-02 8.987799008120435e-01 5.206599198276582e-01
            7.887165251015218e-01 9.592101698204790e-01 5.518095415203144e-01
            2.923023686527778e-01 8.445258972623941e-01 4.816237826105740e-01
            8.394376217243196e-01 3.809652690118381e-01 1.424110712395776e-01
            4.846342847946065e-01 8.960185207187887e-02 3.418619965784917e-01
            2.334785238375495e-01 7.327374078504292e-01 8.649243853491267e-01
            2.205352977884234e-01 6.825735940187665e-01 2.467524046762309e-01
            3.142851290940286e-01 4.351688285299408e-01 1.864740525894623e-01
            9.954273646365726e-01 1.432037715827612e-01 5.262955717159060e-01
            1.891316665873515e-01 3.071053320111327e-01 3.826073822533815e-01
            8.504909683614519e-01 6.803165945287264e-01 8.443289313783434e-01
            8.439047969320805e-01 5.039450455802463e-01 7.090100975839108e-01
            5.283165093245552e-01 7.346571440404910e-01 6.851427434058195e-01
            6.648713312821694e-01 6.480202335624049e-01 3.098340674627222e-01
            3.696163705442677e-01 6.783099141041915e-02 1.444455982426607e-02
            3.174375253775342e-01 9.946974392286094e-01 3.235770414225648e-01
            8.103239146721836e-01 1.811526620489341e-01 7.495520393056689e-01
            1.026232986738559e-01 4.604561345245859e-01 2.231975304627742e-01
            2.127629226746142e-01 5.607414486991632e-01 9.553640662856649e-01
            4.464221878981485e-01 6.561928271414434e-01 1.556289076606221e-02
            3.865511714825203e-01 8.122197154981052e-01 8.123942529260979e-01
            7.746478460871526e-01 7.763636684981615e-01 1.971262118494637e-01
            4.742311326067040e-01 2.965723192828372e-01 2.539214418365823e-01
            7.967912193442113e-01 1.802371328932153e-01 4.496548804382030e-01
            2.485947332486738e-02 8.419145008426521e-01 5.826603970468174e-01
            9.991400610693252e-01 6.544560771481346e-01 5.937059783087021e-01
            2.830534311451685e-01 6.079378808146357e-01 9.670025751544015e-01
            7.678249183946985e-01 4.774558409143720e-01 4.760521844080304e-01
            3.506371321233281e-02 6.115789035631233e-01 6.218175829763063e-01
            8.703594623155821e-01 3.782814048194011e-01 6.172120561564026e-01
            7.429449593511781e-01 8.663442057643033e-01 6.175784189776810e-01
            1.092930448682813e-01 4.802637566363623e-01 3.740494190473312e-01
            8.666167360444277e-01 4.610176534875214e-01 6.945724673749079e-01
            8.215255701244644e-01 3.016271156862274e-01 3.046774258180375e-01
            3.791774365679346e-03 5.546959474092793e-01 3.787943611265603e-01
            6.602632833654598e-01 5.759602723163398e-01 7.566546364987188e-01
            6.323649306580985e-01 1.441944639670529e-01 5.070253148574402e-01
            2.104814884097670e-01 7.212245969560701e-01 6.150490190713972e-01
            4.349466493336438e-01 2.000868787370786e-01 4.020826139637483e-01
            9.331094011171007e-01 2.795344348140895e-01 1.979948684184137e-01
            4.572827278747372e-01 6.986436642689750e-01 5.042247473924152e-01
            1.100286033637933e-01 1.906042478977387e-01 4.701139032209813e-01
            6.798913224075467e-01 4.219296923406707e-01 1.392397488596244e-01
            9.717322457613400e-01 2.076886363023375e-01 4.677804560962281e-01
            4.353773148508129e-01 1.386866282469295e-01 4.799818822685469e-01
            2.379589231318251e-01 8.137607812568703e-01 4.216016518106014e-02
            8.874207766059145e-01 3.284571473601868e-01 9.128779850584334e-02
            9.527711918452563e-01 4.805853317780656e-01 3.183965567295159e-01
            3.916549252268440e-01 6.699177086993990e-02 9.953658473632414e-01
            7.633280411754474e-01 8.204607475443504e-02 6.183741810662485e-01
            1.000527505981699e-01 3.174592244649934e-01 7.557754590909408e-01
            2.238621216085058e-01 3.145782384893810e-01 2.761218757599072e-01
            1.024774992629658e-01 3.380945875697038e-01 5.339715651380794e-01
            2.142099622889593e-02 6.540301366787659e-01 2.258581139947069e-02
            5.169282726791766e-01 4.202204323462601e-01 8.234979057465797e-01
            2.586861670936028e-01 8.320841470115895e-01 2.409693961047334e-01
            3.390302696070702e-02 8.742144834041427e-01 6.521054956857485e-01
            7.496818812261054e-01 5.496928578214269e-01 6.289116344171090e-01
            5.417122635334792e-01 3.380142789846545e-01 1.596561782544259e-01
            8.819482287288631e-01 5.419828098092847e-01 7.571034866051984e-01
            5.613907344033219e-01 2.507236237051315e-01 3.318369889693418e-01
            9.686901623591658e-01 3.681723336565869e-01 2.362803660310983e-01
            7.524751522959622e-01 8.497167766802085e-01 7.107767835752652e-01
            9.492829448715467e-01 8.652630018575396e-01 3.272748519796056e-01
            7.388139122696784e-01 7.411257178158415e-02 5.151409914486603e-01
            2.635250141121041e-01 9.449848289383427e-01 4.321100338515496e-01
            8.748457093283835e-02 1.531741304185485e-02 7.987924639822017e-01
            6.858846169729618e-01 7.927183346497735e-02 5.428664322869750e-01
            8.104648698294635e-01 7.538444689770972e-01 9.178357279881664e-01
            1.508072859583468e-01 3.728655034157864e-01 1.392744910005675e-01
            3.138136268012272e-01 2.094255103101885e-01 9.015069712177659e-01
            4.870284202822635e-01 8.362550723152278e-01 2.762815826771385e-01
            7.036310243232798e-01 8.580924067181045e-01 9.436854303624758e-01
            5.960852029618140e-01 7.727236363010035e-01 6.504108258285968e-02
            5.435115517540157e-02 3.480915769014549e-01 5.509876998454153e-01
            9.711614415409076e-01 9.676103935817125e-01 1.211912112502667e-02
            5.743721802843332e-01 3.250958498358480e-01 1.400148357589668e-01
            5.744430153193537e-01 5.641085741368064e-01 8.560402933234028e-01
            5.771541053555250e-01 2.213889248612743e-01 6.386696038669383e-01
            6.009199560089991e-01 5.870843466572438e-01 8.298836396593373e-01
            1.360163466016264e-01 3.167072840896493e-01 6.990623273340048e-01
            3.493299511449055e-02 7.049518769073405e-01 1.397856742277211e-01
            4.138368189052062e-01 7.176866669684911e-02 6.689369778808085e-01
            2.729493242564299e-01 9.103438573321987e-02 4.666776885056663e-01
            6.896484256476165e-01 1.847378049574832e-01 9.259487076365601e-01
            9.855907923906249e-01 1.925676936090019e-01 8.851382663623240e-01
            2.469153881801451e-01 1.564730728554724e-01 7.069609378974520e-01
            7.663483509827098e-01 2.723621153607964e-01 6.434024616589595e-01
            6.666598346685667e-01 1.827186264538636e-01 2.565425050086261e-01
            5.363164379377821e-01 3.530099378451212e-02 7.035474823697898e-01
            9.005246175085582e-01 7.096017806918975e-02 7.356166435767927e-01
            6.669205235628356e-01 2.244717325446044e-02 7.637923779097499e-01
            5.433847783634581e-01 6.132746315024059e-01 3.947505820585510e-01
            8.927181748011814e-02 2.709410332045040e-01 6.227516046298487e-02
            2.778825256731465e-02 4.696738110753019e-01 7.002586761879319e-01
            9.597927814794447e-01 1.757209246802928e-01 3.956146530656729e-01
            7.322431326499556e-01 8.110932213016965e-01 8.355857555220528e-01
            8.906599287283194e-01 1.933158103752142e-02 2.442230670682800e-01
            8.952530739813831e-01 5.132458435963937e-01 9.612143970031523e-01
            7.943243168265935e-01 9.177456920468011e-01 9.067250774398845e-01
            5.352176744349137e-01 7.482181538304939e-01 3.187265466545215e-01
            3.424007384189467e-01 7.796458549924021e-01 7.672447504611896e-01
            2.807844898412100e-01 3.654713294806605e-01 7.554094546289455e-01
            4.744614549767701e-01 9.631554930981563e-01 5.965433316077267e-01
            8.579637834635889e-03 3.549685478995143e-01 1.427278321854423e-01
            4.248002446727603e-01 8.735699087963896e-01 4.690651194951290e-02
            7.179867474748483e-01 7.348597347267654e-01 8.629947124695373e-01
            4.105359999671968e-01 4.223758175790782e-01 5.924630667145112e-01
            6.045711571835829e-01 5.142482979352739e-01 7.965036990442520e-02
            5.908823645241433e-01 3.358344007026538e-01 2.471112074198448e-01
            9.786060888067165e-01 1.145582301050090e-01 8.469794501660884e-01
            2.759770889324437e-01 7.867557059537083e-01 7.229032433167337e-02
            4.278613978760315e-02 4.753358779104078e-01 8.145752959077540e-01
            8.706533546511643e-01 1.436295280347113e-01 8.843135487104280e-01
            8.878443705041896e-01 4.864560409936232e-01 3.702446220067049e-01
            3.719505800788604e-01 9.069388761651226e-01 8.003447702912697e-01
            1.488739353848852e-01 4.395818304174273e-01 8.106747146687348e-01
            1.693152413252820e-01 6.575470671604471e-01 1.785504435837865e-01
            1.628228890155541e-01 4.605651333802623e-01 6.942936861934366e-01
            5.952206248146819e-01 1.159123772831533e-01 6.609771766172823e-01
            8.514655381871059e-01 3.491162256487661e-01 9.900625435987682e-01
            4.411283131720589e-02 2.978382858106998e-01 3.835147521329627e-01
            5.901110103322387e-01 5.059037324002813e-03 7.305051814818380e-01
            2.477925249665496e-01 1.665762063617485e-01 7.619928942833754e-02
            2.184888674242261e-01 8.745910763757698e-01 7.462615740441292e-02
            7.284449937286197e-01 7.350922459607900e-01 5.348849796022819e-01
            1.830815016582249e-01 6.989792108119663e-02 3.848443335028362e-01
            1.069376438197015e-01 1.492263889465356e-01 3.752260180303243e-01
            4.531075565705437e-01 9.396158033595153e-01 9.729196742043152e-01
            1.937152617774174e-01 8.809344570106784e-01 7.130749127463079e-01
            3.719590499146106e-01 3.246725692692479e-01 3.458911085915399e-01
            1.837605679561112e-01 4.632282047711485e-01 2.437550372010706e-01
            8.867692446644989e-01 7.178294379473661e-01 7.556579255204443e-01
            2.134866119721702e-01 4.284490246738184e-01 8.848847007113212e-01
            5.427747918961460e-01 5.359415547591522e-01 3.465730287881800e-01
            1.652750718712289e-01 9.502837417339915e-01 7.341437406457885e-02
            4.482726416889655e-01 4.980410126779095e-01 5.527935759757688e-01
            8.226081560074687e-01 8.175093958568160e-01 4.166894181073011e-01
            4.725411346371838e-01 7.357272930590890e-01 2.244737724695336e-01
            8.795064946093361e-01 1.581991780443959e-01 4.782280912420008e-01
            6.564112726430429e-01 7.916248798564841e-01 4.993051849896836e-01
            5.479092423818547e-01 1.368335678399923e-02 6.029136853771837e-01
            5.185608344228649e-01 7.234894126151813e-01 2.858026211855519e-01
            4.253576300612313e-01 4.924428442397283e-01 5.395706026641341e-01
            8.121188172426806e-01 3.374043796001711e-01 5.195122910393892e-01
            6.420299388128383e-01 7.548169927084709e-01 4.781912695979231e-01
            1.842726229091473e-01 7.064911374633832e-01 5.809584709966044e-01
            9.977038385029215e-01 4.425615921879329e-01 3.908356321060984e-01
            4.883876151919846e-01 8.215931159897113e-01 9.205608556235451e-01
            5.092354685943664e-01 8.160267487608520e-01 6.763904281396230e-01
            7.641759992232936e-01 8.161730117135719e-01 2.914656229354518e-01
            6.745894957706815e-01 5.010604959659127e-01 9.636297267937333e-01
            7.252916466095685e-01 1.637207116450736e-01 1.638846617985701e-02
            2.222603380251766e-01 5.840228108371669e-01 2.220669994947079e-01
            3.830849388007008e-01 9.091202425875744e-01 8.287269124358693e-01
            3.461632147477373e-01 2.896303253705851e-01 7.788367716944967e-01
            ];
        
        
        alphashape_radius=[
            5.000000000000000e-01
            ];
        
        
        
    case -7
        
        polyhedron_type='Concave Polyhedron: one hole';
        
        vertices=[
            8.682329144425175e-01 3.266410024432838e-01 5.750103154741084e-01
            8.186683735954834e-01 3.251516940049070e-01 7.411378921030854e-01
            5.900841322901613e-01 4.001734208480078e-01 7.536707305663618e-01
            2.287841210027897e-03 7.635061191199046e-01 5.550881101582001e-01
            9.445320434483669e-01 1.640677875170120e-01 8.009122590275306e-01
            7.924769014400519e-01 6.882903484206540e-01 6.703459822631603e-01
            3.993284292658533e-01 6.433166187171735e-01 5.351074809664100e-01
            8.734494981429497e-01 8.213535257649524e-01 6.859261164376368e-01
            6.743963258964303e-01 3.785469303302259e-01 9.551838395988530e-01
            2.423656059243389e-01 1.119006607697201e-01 6.277522906099283e-01
            1.684447239530630e-01 8.594832743733937e-01 8.191319965204129e-01
            8.611248013709782e-01 7.131831063374162e-01 1.889653326543662e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
        
    case -6
        
        polyhedron_type='Concave Polyhedron; 2 disconnected polyhedra';
        
        vertices=[
            7.043122976073990e-01 7.470176965745213e-01 7.730979700393871e-01
            4.970221828134525e-01 3.099280908647454e-01 9.319814294926577e-01
            6.095968655514951e-01 9.165041110106615e-01 1.768844290079774e-01
            7.782840158881210e-02 7.822659721629751e-01 2.042176149515484e-02
            2.414326962793247e-01 1.010467745665099e-01 8.336389206206382e-01
            2.578415583689841e-02 2.811045519926746e-01 4.731935927192107e-01
            8.464769359226610e-01 9.075052740697294e-01 1.086462702273944e-01
            1.486279613221164e-02 1.110871585270446e-01 8.905437818186224e-02
            5.504158884128526e-01 7.898195137681069e-01 1.151295876565336e-01
            1.196711575040840e-01 1.022795552424890e-01 1.679728584195933e-01
            5.879491266122041e-01 7.679438903195279e-01 9.725103383139468e-01
            6.222109441717147e-01 6.620554240394096e-01 2.860036350944583e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
    case -5
        
        polyhedron_type='Concave Polyhedron; 3 polyhedra with one vertex or side in common (resp.)';
        
        vertices=[
            8.001115866786623e-02 6.376410230785723e-03 7.051649108289371e-02
            7.935380917331858e-02 9.199013631348676e-01 6.002365591523540e-02
            7.237624919866215e-02 1.796224130559454e-02 6.615771015291297e-01
            9.003410060810889e-01 2.941647987706497e-02 3.440577694777694e-01
            9.693082988664113e-01 7.114519163682994e-01 1.285133482338328e-01
            3.913069588115609e-01 5.384141865401406e-01 6.371356638478847e-01
            3.135851048672553e-01 5.245885805189759e-01 7.464881699547887e-01
            5.533135799374451e-01 5.022399304161738e-01 8.052672976354331e-01
            7.920480624973907e-01 6.567249822527477e-02 2.690271315860681e-01
            7.983140830802098e-01 2.300406237580992e-01 4.339797789833593e-01
            8.632585405599013e-01 1.170403198125443e-01 4.018064872641355e-01
            7.979622044473094e-01 9.897710614201980e-01 6.081423394349520e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
        
    case -4
        
        polyhedron_type='Concave Polyhedron; 2 polyhedra with two vertices in common and hole';
        
        vertices=[
            7.019456892102048e-01 6.559481695335554e-01 8.318498806648059e-01
            3.928739144202648e-01 8.565376013216209e-01 3.944709820411555e-01
            5.718686630258228e-01 2.709404658120863e-01 9.663634688171530e-02
            2.497737351837366e-01 8.959477796220956e-01 2.653916148253564e-01
            7.273467331114006e-01 2.087903958723814e-01 7.919450125521071e-01
            3.433291453593958e-01 9.641681802460583e-01 9.368607426003079e-01
            7.863467486325948e-02 4.143600880391674e-01 5.682585883425429e-01
            3.022309282744359e-01 2.295652138635622e-01 9.380215572322496e-01
            3.004024819239974e-01 8.183701989647573e-01 5.972483945500038e-01
            2.322119759663466e-02 3.604434932620764e-01 9.879519676560238e-01
            1.761142516577059e-01 5.015240681494324e-01 7.623485025252664e-01
            8.609187966110264e-01 7.476741184362131e-01 1.856158686068560e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
        
    case -3
        
        polyhedron_type='Concave Polyhedron; 3 polyhedra with one vertex in common';
        
        vertices=[
            5.929069367828630e-01 1.176061171303259e-01 4.559313543318710e-02
            6.576080119544996e-01 2.613032698813186e-01 6.554734759809633e-01
            3.489626928986161e-01 8.933580019130812e-01 7.810422089046415e-01
            5.369344743552795e-02 9.720611376789213e-01 8.629004252215658e-01
            4.716803107485978e-02 9.940738862116841e-01 7.853528718331544e-01
            5.585780916850613e-02 3.467862632547470e-01 4.334923027234355e-01
            7.151199562065402e-02 2.923021355683054e-01 3.478212469654292e-01
            1.567554495465118e-02 1.872474071594036e-01 5.538855025922651e-02
            9.736790072546389e-01 7.431483363186397e-01 8.194414425887224e-01
            7.347967209312096e-01 7.180187226846149e-01 1.060403382092259e-01
            7.523856147363535e-01 1.888947255123814e-01 6.074046531967858e-03
            3.331360813113482e-02 4.339378046020201e-01 2.750423525483723e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
        
    case -2
        
        polyhedron_type='Concave Polyhedron; 2 polyhedra with one vertex in common';
        
        vertices=[
            8.605621147246640e-01 4.478896435877100e-01 9.503684015148999e-01
            9.390935153524338e-01 6.512486632169324e-01 6.942862475655555e-01
            6.801939644480517e-01 1.695018706149451e-01 2.068065270209161e-01
            9.174237484474340e-01 5.314492413595803e-01 5.547623259567388e-01
            2.566916709011579e-01 6.338008581939598e-01 8.792784845119187e-01
            8.856181344290193e-01 1.409600418269241e-02 5.578576907197960e-01
            9.200427901181818e-01 4.703713352310025e-01 7.523340601514115e-01
            3.000634400366826e-01 8.863259711681059e-01 8.949007922295080e-01
            7.339094669379265e-02 1.140284356780709e-01 8.418441730648902e-01
            7.673987349601374e-01 4.425411516501793e-01 1.308566826459481e-01
            8.495221586501400e-02 6.595478881615526e-01 1.891539338519886e-01
            7.287640795926718e-01 2.947729702676847e-01 1.536398308925564e-01
            ];
        
        
        alphashape_radius=[
            4.000000000000000e-01
            ];
        
    case -1
        polyhedron_type='Polyhedron';
        numvertices=10;
        vertices=rand(numvertices,3);
        alphashape_radius=0.5;
        
        FID=fopen('vertices.txt','w');
        fprintf(FID,'\n vertices=[');
        fprintf(FID,'\n %1.15e %1.15e %1.15e',vertices');
        fprintf(FID,'\n ]; \n \n');
        
        fprintf(FID,'\n alphashape_radius=[');
        fprintf(FID,'\n %1.15e',alphashape_radius);
        fprintf(FID,'\n ]; \n \n');
        
        fprintf(FID,'\n \n %% ------------------------------------------------------  \n \n');
        
        s=fclose(FID);
        
        
        
    case 0
        polyhedron_type='Polyhedron';
        numvertices=15;
        vertices=rand(numvertices,3);
    case 1 % Unit simplex
        polyhedron_type='Unit simplex';
        A=[0 0 0]; B=[1 0 0]; C=[0 1 0]; D=[0 0 1];
        vertices=[A; B; C; D];
    case 2 % Unit tetrahedron
        polyhedron_type='Unit tetrahedron';
        A=[1 -1/sqrt(3) -1/sqrt(6)];
        B=[-1 -1/sqrt(3) -1/sqrt(6)];
        C=[0 2/sqrt(3) -1/sqrt(6)];
        D=[0 0 3/sqrt(6)];
        vertices=[A; B; C; D];
    case 3 % Unit tetrahedron via cells
        polyhedron_type='Unit tetrahedron (via cells)';
        A=[1 -1/sqrt(3) -1/sqrt(6)];
        B=[-1 -1/sqrt(3) -1/sqrt(6)];
        C=[0 2/sqrt(3) -1/sqrt(6)];
        D=[0 0 3/sqrt(6)];
        vertices=[A; B; C; D];
    case 4
        polyhedron_type='Tetrahedron (shifted)';
        c=[1 1 1];
        % A=[0 0 0]; B=[1 0 0]; C=[0 1 0]; D=[0 0 1];
        A=[-1 -1/sqrt(3) -1/sqrt(6)]; B=[0 2/sqrt(3) -1/sqrt(6)];
        C=[1 -1/sqrt(3) -1/sqrt(6)]; D=[0 0 3/sqrt(6)];
        A=c+A; B=c+B; C=c+C; D=c+D;
        vertices=[A; B; C; D];
    case 5 % Unit cube
        polyhedron_type='Unit cube (via matrix)';
        A=[0 0 0]; B=[1 0 0]; C=[1 1 0]; D=[0 1 0];
        A1=[0 0 1]; B1=[1 0 1]; C1=[1 1 1]; D1=[0 1 1];
        vertices=[A; B; C; D; A1; B1; C1; D1];
    case 6 % Unit cube via cells
        polyhedron_type='Unit cube (via cells)';
        A=[0 0 0]; B=[1 0 0]; C=[1 1 0]; D=[0 1 0];
        A1=[0 0 1]; B1=[1 0 1]; C1=[1 1 1]; D1=[0 1 1];
        vertices=[A; B; C; D; A1; B1; C1; D1];
    case 7 % Shifted Unit cube
        polyhedron_type='Shifted unit cube';
        A=[0 0 0]; B=[1 0 0]; C=[1 1 0]; D=[0 1 0];
        A1=[0 0 1]; B1=[1 0 1]; C1=[1 1 1]; D1=[0 1 1];
        vertices=[A; B; C; D; A1; B1; C1; D1]-0.5;
    case 8
        polyhedron_type='Tetrahedron';
        vertices=[2.991748564344140e-01     3.012055865134726e-01     5.821813351685147e-01
            5.311265268153221e-01     9.505815914084504e-01     1.530686974430926e-01
            1.463175172423004e-03     4.606481932062152e-01     7.309430616949053e-02
            8.837888094022676e-01     2.876475110417300e-01     5.805647768308455e-01];
    case 9
        polyhedron_type='Tetrahedron';
        vertices=[5.311265268153221e-01     9.505815914084504e-01     1.530686974430926e-01
            4.043838573077683e-01     8.462663857638142e-02     2.870147523971093e-01
            1.463175172423004e-03     4.606481932062152e-01     7.309430616949053e-02
            8.837888094022676e-01     2.876475110417300e-01     5.805647768308455e-01];
    case 10
        polyhedron_type='Spherical polyhedron (4)';
        [x1, y1, z1] = sphere(4);
        P=[x1(:) y1(:) z1(:)];
        vertices = unique(P,'rows');
        alphashape_radius=realmax;
    case 11
        polyhedron_type='Spherical polyhedron (10)';
        [x1, y1, z1] = sphere(10);
        P=[x1(:) y1(:) z1(:)];
        vertices = unique(P,'rows');
        alphashape_radius=realmax;
    case 12
        polyhedron_type='Spherical polyhedron (20)';
        [x1, y1, z1] = sphere(20);
        P=[x1(:) y1(:) z1(:)];
        vertices = unique(P,'rows');
        alphashape_radius=realmax;
    case 13
        polyhedron_type='Spherical polyhedron (40)';
        [x1, y1, z1] = sphere(40);
        P=[x1(:) y1(:) z1(:)];
        vertices = unique(P,'rows');
        alphashape_radius=realmax;
end

c=[1 1 1];
% c=[0 0 0];
c=repmat(c,size(vertices,1),1);
vertices=c+vertices;










function moments=evaluate_monomial_moments(XYZW,ade)

%--------------------------------------------------------------------------
% Object:
% Moments of the rule with nodes XYZW(:,1:3) and weights XYZW(:,4),
% relative to monomial basis of total degree "ade" (lexicographic order).
%--------------------------------------------------------------------------
% Input:
% ade: polynomial degree.
%--------------------------------------------------------------------------
% Output:
% moments: approximation of the moments of the Vandermonde matrix
%    of total degree "ade",
%--------------------------------------------------------------------------
% Authors:
% Alvise Sommariva and Marco Vianello
% University of Padova, September 7, 2020
%--------------------------------------------------------------------------

% ............................ troubleshooting ............................

X=XYZW(:,1); Y=XYZW(:,2); Z=XYZW(:,3); W=XYZW(:,4);

% .......................... moments computation ..........................

nodes=XYZW(:,1:3);


moments=[];
for ii=0:ade
    for jj=0:(ade-ii)
        for kk=0: ade-jj-ii % computing contribution
            FXYZC=(X.^ii).*(Y.^jj).*(Z.^kk);
            moments_L=W'*FXYZC;
            moments=[moments; moments_L];
        end
    end
end










