function [result,m] = validmesh(m,verbose)
%RESULT = VALIDMESH(M,VERBOSE)  Determines whether the mesh is valid.
%  All the connectivity information is checked for consistency.
%  Triangles with very small angles are detected.
%  The morphogen array must be the right size.
%  If VERBOSE is true (default false), messages will be printed to the
%  command window in case of flaws.
%  RESULT will be FALSE if the mesh has fatal errors in it.

if isempty(m), result = 0; return; end

if nargin < 2
    verbose = false;
end
ERRORS = false;
if ERRORS
    complainer = @error;
else
    complainer = @warning;
end
if verbose
    whiner = @warning;
else
    whiner = @donothing;
end
    
result = 1;
if isfield( m, 'globalProps' ) ...
        && isfield( m.globalProps, 'validateMesh' ) ...
        && ~m.globalProps.validateMesh
    return;
end

setGlobals();
global gDYNAMICFIELDS gSTATICFIELDS gHYBRIDFIELDS gUNSAVEDFIELDS gTRANSIENTFIELDS
[ok,missingfields,extrafields] = checkFields( m, ...
    { gDYNAMICFIELDS{:} gSTATICFIELDS{:} gHYBRIDFIELDS{:} gUNSAVEDFIELDS{:} }, ...
    gTRANSIENTFIELDS);
if ~isempty(missingfields)
    complainer( 'Mesh has missing fields:' );
    for i=1:length(missingfields)
        fprintf( 1, '\n    %s', missingfields{i} );
    end
    fprintf( 1, '\n' );
end
if ~isempty(extrafields)
    complainer( 'Mesh has extra fields:' );
    for i=1:length(extrafields)
        fprintf( 1, '\n    %s', extrafields{i} );
    end
    fprintf( 1, '\n' );
end

np = size(m.nodes,1);
ne = size(m.edgeends,1);
nf = size(m.tricellvxs,1);

if ~checknumel( m, 'celldata', nf, complainer );
    result = 0;
end

% Check that every triangle has three distinct vertexes.
for fi=1:nf
    vxs = m.tricellvxs(fi,:);
    if (vxs(1)==vxs(2)) || (vxs(3)==vxs(1)) || (vxs(2)==vxs(3))
        result = 0;
        complainer( 'Triangle %d fails to have distinct vertexes: [%d,%d,%d].\n', ...
            fi, vxs(1), vxs(2), vxs(3) );
    end
end

% Check that every triangle has three distinct edges.
for fi=1:nf
    es = m.celledges(fi,:);
    if (es(1)==es(2)) || (es(3)==es(1)) || (es(2)==es(3))
        result = 0;
        complainer( 'Triangle %d fails to have distinct edges: [%d,%d,%d].\n', ...
            fi, es(1), es(2), es(3) );
    end
end

% Check that for every connection from a cell to an edge in celledges,
% the edge index is valid, and there is a connection from that edge to that
% cell in edgecells.
% if any( all( m.celledges ~= repmat((1:nf)',1,3), 2 ), 1 )
    for fi=1:nf
        for ei=m.celledges(fi,:)
            if all(m.edgecells(ei,:) ~= fi)
                result = 0;
                complainer( 'validmesh:celledge2', ...
                    'Cell %d has edge %d but edge is connected to cells %d and %d.', ...
                    fi, ei, m.edgecells(ei,1), m.edgecells(ei,2) );
            end
        end
    end
% end
if any(any((m.celledges<1) | (m.celledges > ne)))
    for fi=1:nf
        eis = m.celledges(fi,:);
        for ei=eis
            if (ei < 1) || (ei > ne)
                result = 0;
                complainer( 'validmesh:celledge1', ...
                    'Cell %d has edge %d, which is outside the valid range %d:%d.', ...
                    fi, ei, 1, ne );
            end
        end
    end
end

% Check that the ends of each edge and the cells on either side of the
% edge, as given by edgeends and edgecells, are valid indexes.
% Check that the ends of each edge are vertexes of the cells of that edge.
% Check that the edge ends are distinct.
% Check that every edge has a cell on at least one side.
% Check that if there is a cell on both sides of the edge, the two cells
% are distinct.
if any(any(m.edgeends<1))
    result = 0;
    for ei=1:ne
        p1 = m.edgeends(ei,1);
        if p1 < 1
            complainer( 'validmesh:badpointindex1', ...
                'Edge %d endpoint 1 is %d, should be positive.', ...
                ei, p1 );
        end
        p2 = m.edgeends(ei,2);
        if p2 < 1
            complainer( 'validmesh:badpointindex2', ...
                'Edge %d endpoint 2 is %d, should be positive.', ...
                ei, p2 );
        end
    end
end
sameends = m.edgeends(:,1)==m.edgeends(:,2);
if any(sameends)
    result = 0;
    for ei=1:ne
        if sameends(ei)
            complainer( 'validmesh:badpointindex3', ...
                'Edge %d endpoints should be distinct, are %d %d', ...
                ei, m.edgeends(ei,:) );
        end
    end
end
if any(any(m.edgeends>np))
    result = 0;
    for ei=1:ne
        p1 = m.edgeends(ei,1);
        if p1 > np
            complainer( 'validmesh:badpointindex1', ...
                'Edge %d endpoint 1 is %d, should be not more than %d.', ...
                ei, p1, np );
        end
        p2 = m.edgeends(ei,2);
        if p2 > np
            complainer( 'validmesh:badpointindex2', ...
                'Edge %d endpoint 2 is %d, should be not more than %d.', ...
                ei, p2, np );
        end
    end
end
badcells = m.edgecells(:,1) < 1;
if any(badcells)
    result = 0;
    for ei=find(badcells)'
        complainer( 'validmesh:badpointindex4', ...
            'Edge %d cell 1 is %d, should be positive.', ...
            ei, m.edgecells(ei,1) );
    end
end
badcells = m.edgecells(:,1) < 0;
if any(badcells)
    for ei=find(badcells)'
        result = 0;
        complainer( 'validmesh:badpointindex4', ...
            'Edge %d cell 2 is %d, should be non-negative.', ...
            ei, m.edgecells(ei,2) );
    end
end
badcells = m.edgecells > nf;
if any(badcells(:))
    result = 0;
    for ei=1:nf
        for j=1:2
            complainer( 'validmesh:badpointindex4', ...
                'Edge %d cell %d is %d, should be no more than %d.', ...
                ei, j, m.edgecells(ei,j), nf );
        end
    end
end
badcells = m.edgecells(:,1)==m.edgecells(:,2);
if any(badcells(:))
    result = 0;
    for ei=find(badcells)'
        complainer( 'validmesh:badpointindex6', ...
            'Edge %d cells on either side should be distinct, are %d %d', ...
            ei, m.edgecells(ei,:) );
    end
end
for ei=1:ne
    f1 = m.edgecells(ei,1);
    f2 = m.edgecells(ei,2);
    if (f1 ~= 0) && ~isedgeof( m, ei, f1 )
        result = 0;
        complainer( 'validmesh:badpointindex7', ...
            'Edge %d (n %d,%d) (f %d,%d) is not an edge or not the right edge of cell %d (p %d,%d,%d, e %d,%d,%d)', ...
            ei, m.edgeends(ei,1), m.edgeends(ei,2), ...
            f1, f2, f1, ...
            m.tricellvxs(f1,1), m.tricellvxs(f1,2), m.tricellvxs(f1,3), ...
            m.celledges(f1,1), m.celledges(f1,2), m.celledges(f1,3) );
    end
    if (f2 ~= 0) && ~isedgeof( m, ei, f2 )
        result = 0;
        complainer( 'validmesh:badpointindex8', ...
            'Edge %d (n %d,%d) (f %d,%d) is not an edge or not the right edge of cell %d (p %d,%d,%d, e %d,%d,%d)', ...
            ei, m.edgeends(ei,1), m.edgeends(ei,2), f1, f2, f2, ...
            m.tricellvxs(f2,1), m.tricellvxs(f2,2), m.tricellvxs(f2,3), ...
            m.celledges(f2,1), m.celledges(f2,2), m.celledges(f2,3) );
    end
end

% Check that two cells do not have the same vertexes.
    errs = checkUniqueRows( m.tricellvxs );
    if ~isempty(errs)
        result = 0;
        complainer( 'validmesh:dupcellvxs', ...
            'Some cells have the same vertexes:' );
        fprintf( 1, ' %d', errs );
        fprintf( 1, '\n' );
    end
% Check that two cells do not have the same edges.
    errs = checkUniqueRows( m.celledges );
    if ~isempty(errs)
        result = 0;
        complainer( 'validmesh:dupcelledges', ...
            'Some cells have the same edges:' );
        fprintf( 1, ' %d', errs );
        fprintf( 1, '\n' );
    end
% Check that two edges do not have the same ends.
    errs = checkUniqueRows( m.edgeends );
    if ~isempty(errs)
        result = 0;
        complainer( 'validmesh:dupedgeends', ...
            'Some edges have the same ends:' );
        fprintf( 1, ' %d', errs );
        fprintf( 1, '\n' );
    end

% Check that the mesh is oriented.
% This means that every edge that belongs to two cells should have its ends
% occur in opposite orders in those cells.
checkOrientation = 1;
if checkOrientation
    for ei=1:ne
        c2 = m.edgecells(ei,2);
        if c2==0, continue; end
        c1 = m.edgecells(ei,1);
        n1 = m.edgeends(ei,1);
        n2 = m.edgeends(ei,2);
        c1n1 = find(m.tricellvxs(c1,:)==n1);
        c1n2 = find(m.tricellvxs(c1,:)==n2);
        c2n1 = find(m.tricellvxs(c2,:)==n1);
        c2n2 = find(m.tricellvxs(c2,:)==n2);
        fwd1 = mod(c1n2-c1n1,3)==1;
        fwd2 = mod(c2n2-c2n1,3)==1;
        if fwd1==fwd2
            % Orientation error
            result = 0;
            complainer( 'validmesh:orientation', ...
                'Cells %d [%d %d %d] and %d [%d %d %d] contain edge %d with nodes %d and %d in the same order.\n', ...
                c1, m.tricellvxs(c1,1), m.tricellvxs(c1,2), m.tricellvxs(c1,3), ...
                c2, m.tricellvxs(c2,1), m.tricellvxs(c2,2), m.tricellvxs(c2,3), ...
                ei, n1, n2 );
        end
    end
end

% Check that every node belongs to some cell.
allCellNodes = unique(reshape(m.tricellvxs,1,[]));
checkallints( 'validmesh tricellvxs', allCellNodes, np );
% Check that every edge belongs to some cell.
allCellEdges = unique(reshape(m.celledges,1,[]));
checkallints( 'validmesh celledges', allCellEdges, ne );
% Check that every node belongs to some edge.
allEdgeEnds = unique(reshape(m.edgeends,1,[]));
checkallints( 'validmesh edgeends', allEdgeEnds, np );
% Check that every cell belongs to some edge.
allEdgeCells = unique(reshape(m.edgecells,1,[]));
if (~isempty(allEdgeCells)) && (allEdgeCells(1)==0)
    allEdgeCells = allEdgeCells(2:end);
end
checkallints( 'validmesh edgecells', allEdgeCells, nf );

% Check that nodecelledges is valid.
badnodecelledges = false;
if ~isfield( m, 'nodecelledges' )
    badnodecelledges = true;
    whiner( 'validmesh:nodecelledges', ...
        'nodecelledges field is missing.\n' );
else
    if length(m.nodecelledges) ~= size(m.nodes,1)
        badnodecelledges = true;
        whiner( 'validmesh:nodecelledges', ...
            'nodecelledges has the wrong length: %d found, but %d nodes in the mesh.\n', ...
            length(m.nodecelledges), size(m.nodes,1) );
    else
        for vi = 1:length(m.nodecelledges)
            nce = m.nodecelledges{vi};
            if isempty(nce)
                warning( 'validmesh:chains1', ...
                    'Node %d has an empty array of neighbours.', ...
                    vi );
                ok = false;
            else
                nbedges = nce(1,:);
                nbcells = nce(2,:);
                % Check that vi is an end of every edge in nbedges
                badnbedges = find( ~any( m.edgeends( nbedges, : )==vi, 2 ) );
                if any( badnbedges )
                    badnodecelledges = true;
                    whiner( 'validmesh:nodecelledges', ...
                        'node %d is not a neighbour of edges %s which are in its neighbour list %s.\n', ...
                        vi, ...
                        nums2string( nbedges(badnbedges), '%d' ), ...
                        nums2string( nbedges, '%d' ) );
                end
                nznbcells = find( nbcells ~= 0 );
                badnbcells = find( ~any( m.tricellvxs( nbcells(nznbcells), : )==vi, 2 ) );
                if any(badnbcells)
                    badnodecelledges = true;
                    whiner( 'validmesh:nodecelledges', ...
                        'node %d is not a vertex of cells %s which are in its neighbour list %s.\n', ...
                        vi, ...
                        nums2string( nbcells(nznbcells(badnbcells)), '%d' ), ...
                        nums2string( nbcells, '%d' ) );
                end
            end
        end
    end
end
if badnodecelledges && result
    m = makeVertexConnections( m );
    complainer( 'validmesh:nodecelledges', ...
        'The nodecelledges structure has been repaired.' );
end

% The angles of each FEM element should not be very small or very close to
% 180 degrees.
ANGLE_THRESHOLD = m.globalProps.mincellangle/4;  % 4 is a fudge factor.
alltriangles = permute( reshape( m.nodes( m.tricellvxs', : ), 3, [], 3 ), [1 3 2] );
allAngles = trianglesAngles( alltriangles );
badAngles = allAngles < ANGLE_THRESHOLD;
if any(badAngles(:))
    complainer( 'validmesh:poorcells', ...
        'There are %d very small angles (below %f radians), smallest is %.3g radians.', ...
        sum(badAngles(:)), ANGLE_THRESHOLD, min(abs(allAngles(:))) );
    if false
        for ci=1:nf
            if any( badAngles(:,ci) )
                fprintf( 1, 'Cell %d has angles [%f %f %f], threshold %f.\n', ...
                    ci, allAngles(:,ci), ANGLE_THRESHOLD );
                fprintf( 1, '    %8f   %8f   %8f\n', alltriangles(1,:,ci) );
                fprintf( 1, '    %8f   %8f   %8f\n', alltriangles(2,:,ci) );
                fprintf( 1, '    %8f   %8f   %8f\n', alltriangles(3,:,ci) );
            end
        end
    end
end

% The morphogens must be defined per vertex.
if isfield(m,'morphogens') && (size(m.morphogens,1) ~= size(m.nodes,1))
    result = 0;
    complainer( 'validmesh:badgrowth', ...
        'Morphogens are defined for %d nodes, but the mesh has %d nodes.\n', ...
        size(m.morphogens,1), size(m.nodes,1) );
    delta = size(m.morphogens,1) - size(m.nodes,1);
    if delta > 0
        m.morphogens = m.morphogens( 1:size(m.nodes,1), size( m.morphogens, 2 ) );
    else
        m.morphogens = [ m.morphogens ; zeros( delta, size( m.morphogens, 2 ) ) ];
    end
end

% Check that cellstiffness is the right size.
cssize = size(m.cellstiffness);
if ~samesize(cssize, [6 6 nf])
    complainer( 'validmesh:badcellstiffness', ...
        'cellstiffness has size [%d %d %d], expected [6 6 %d]\n', ...
        cssize(1), cssize(2), cssize(3), nf );
    if all(cssize([1 2])==[6 6])
        if cssize(3) > nf
            m.cellstiffness = m.cellstiffness(:,:,1:nf);
        else
            m.cellstiffness(:,:,(end+1):nf) = ...
                repmat( sum( m.cellstiffness, 3 )/cssize(3), [1, 1, nf-cssize(3)] );
        end
    else
        ok = false;
    end
end
    
% Check that cellbulkmodulus is the right size.
cbmsize = length(m.cellbulkmodulus);
if cbmsize ~= nf
    complainer( 'validmesh:badcellbulkmodulus', ...
        'cellbulkmodulus has size %d, expected %d\n', ...
        cbmsize, nf );
    if cbmsize > nf
        m.cellbulkmodulus = m.cellbulkmodulus(1:nf);
    else
        m.cellbulkmodulus(:,:,(end+1):nf) = sum( m.cellbulkmodulus )/nf;
    end
end
    
% Check that cellpoisson is the right size.
cpsize = length(m.cellpoisson);
if cpsize ~= nf
    complainer( 'validmesh:badcellpoisson', ...
        'cellpoisson has size %d, expected %d\n', ...
        cpsize, nf );
    if cpsize > nf
        m.cellpoisson = m.cellbulkmodulus(1:nf);
    else
        m.cellpoisson(:,:,(end+1):nf) = sum( m.cellpoisson )/nf;
    end
end
    

% Check that no two vertexes are in the same place.
    if false % No, don't bother.
        if verbose
            mul2 = 1.1234;
            mul3 = pi;
            vxhash = m.nodes(:,1) + mul2*m.nodes(:,2) + mul3*m.nodes(:,3);
            [vxsort,sortperm] = sort(vxhash);
            comparevxs = find( vxsort(1:(end-1)) == vxsort(2:end) );
            for i=1:length(comparevxs)-1
                vx1 = sortperm(comparevxs(i));
                vx2 = sortperm(comparevxs(i+1));
                fprintf( 1, 'Vertexes %d and %d appear to be in the same place:\n', ...
                    vx1, vx2 );
                fprintf( 1, '    [%.3f,%.3f,%.3f]    [%.3f,%.3f,%.3f]\n', ...
                    m.nodes(vx1,:), m.nodes(vx2,:) );
            end
        end
    end
    
% If the prism nodes are valid, check there are the right number.
    if m.globalProps.prismnodesvalid
        if size(m.prismnodes,1) ~= size(m.nodes,1)*2
            result = 0;
            complainer( 'validmesh:badprisms', ...
                'Wrong number of prism nodes: %d, expected twice %d.', ...
                size(m.prismnodes,1), size(m.nodes,1) );
        end
    end
    
    if (~isempty(m.displacements)) && ~checkheight( m, 'displacements', size(m.prismnodes,1), complainer );
        result = 0;
    end
  % if ~checkheight( m, 'effectiveGrowthTensor', size(m.tricellvxs,1), complainer );
  %     result = 0;
  % end
  
    if ~checkheight( m, 'fixedDFmap', size(m.nodes,1), complainer );
        result = 0;
    end
    
% Validate the edge/cell chains.
  % ok = validateChains(m);
  % result = result && ok;
    
% Validate the seams
    if length(m.seams) ~= size(m.edgeends,1)
        whiner( 'validmesh:badseams', ...
            'Seams bitmap has %d elements but there are %d edges.', ...
            length(m.seams), size(m.edgeends,1) );
    end
    
% Validate the second layer.
    if hasSecondLayer( m )
        badvxFEMcell = find( ...
            (m.secondlayer.vxFEMcell > size(m.tricellvxs,1))...
            | (m.secondlayer.vxFEMcell <= 0) );
        if ~isempty(badvxFEMcell)
            result = 0;
            complainer( 'validmesh:badvxFEMcell', ...
                '%d invalid values found in secondlayer.vxFEMcell:', ...
                length(badvxFEMcell) );
            fprintf( 1, ' %d', badvxFEMcell );
            fprintf( 1, '\n' );
        end
    end
        
    
% Check that there is the right amount of other 3D data.
    numcells = size(m.tricellvxs,1);
    numedges = size(m.edgeends,1);
    numnodes = size(m.nodes,1);
    result = checksize( numcells, length(m.cellareas), 'cellareas', complainer ) && result;
    result = checksize( numcells, size(m.unitcellnormals,1), 'unitcellnormals', complainer ) && result;
    result = checksize( numedges, length(m.currentbendangle), 'currentbendangle', complainer ) && result;
    result = checksize( numedges, length(m.initialbendangle), 'initialbendangle', complainer ) && result;
    if ~isempty(m.displacements)
        result = checksize( numnodes*2, size(m.displacements,1), 'displacements', complainer ) && result;
    end

    okpol = true;
    okpol = checksize( numcells, size(m.gradpolgrowth,1), 'gradpolgrowth', complainer ) && okpol;
    okpol = checksize( numcells, size(m.polfreeze,1), 'polfreeze', complainer ) && okpol;
    okpol = checksize( numcells, size(m.polfreezebc,1), 'polfreezebc', complainer ) && okpol;
    okpol = checksize( numcells, size(m.polfrozen,1), 'polfrozen', complainer ) && okpol;
    if ~okpol
        complainer( 'validmesh:badpolgrad', 'Polariser gradient errors fixed.\n' );
        m = calcPolGrad( m );
        result = false;
    end
    
    if ~isempty( m.growthanglepervertex )
        result = checksize( numnodes, size(m.growthanglepervertex,1), 'growthanglepervertex', complainer ) && result;
    end
    if ~isempty( m.growthangleperFE )
        result = checksize( numcells, size(m.growthangleperFE,1), 'growthangleperFE', complainer ) && result;
    end
    if ~isempty( m.decorFEs )
        if length(m.decorFEs) ~= size(m.decorBCs,1)
            complainer( 'validmesh:baddecorpts', ...
                'length(m.decorFEs) = %d but size(m.decorBCs,1) = %d.\n', ...
                length(m.decorFEs), size(m.decorBCs,1) );
            m.decorFEs = [];
            m.decorBCs = [];
        end
        if any(m.decorFEs > numcells)
            complainer( 'validmesh:invalid decor FEs:' );
            m.decorFEs(m.decorFEs > numcells);
        end
        if any(m.decorFEs < 1)
            complainer( 'validmesh:invalid decor FEs:' );
            m.decorFEs(m.decorFEs < 5);
        end
    end

% Check the per-morphogen data.
    numMorphogens = size( m.morphogens, 2 );
    global gPerMgenDefaults
    fns = fieldnames( gPerMgenDefaults );
    for i=1:length(fns)
        fn = fns{i};
        numfound = size(m.(fn),2);
        if numfound > numMorphogens
            complainer( ['validmesh:bad',fn], ...
                '%d extra columns found in %s:', ...
                numfound - numMorphogens, fn );
            x = m.(fn);
            m.(fn) = x(:,1:numMorphogens);
        elseif numfound < numMorphogens
            complainer( ['validmesh:bad',fn], ...
                '%d missing columns found in %s:', ...
                numMorphogens - numfound, fn );
            m.(fn) = [ m.(fn), ...
                       repmat( gPerMgenDefaults.(fn), 1, numMorphogens - numfound ) ];
        end
    end
    
    numIndexedMgens = length(m.mgenIndexToName);
    if numIndexedMgens < numMorphogens
        % Well, fuck.
        complainer( 'validmesh:badmorphogennames', ...
            '%d unnamed morphogens found in mgenIndexToName:', ...
            numMorphogens - numfound );
        for i=(numIndexedMgens+1):numMorphogens
            m.mgenIndexToName{i} = sprintf( 'UNK_%03d', i );
        end
        m.mgenNameToIndex = invertNameArray( m.mgenIndexToName );
    elseif numIndexedMgens > numMorphogens
        complainer( 'validmesh:badmorphogennames', ...
            '%d extra morphogens names in mgenIndexToName:', ...
            numMorphogens - numfound );
        m.mgenIndexToName = m.mgenIndexToName( 1:numMorphogens );
        m.mgenNameToIndex = invertNameArray( m.mgenIndexToName );
    end
    m.mgenNameToIndex = invertNameArray( m.mgenIndexToName );
end

function [ok,v] = rectifysize1( expectedsize, v, defaultval, msg, complainer )
    actualsize = length(v);
    ok = checksize( expectedsize, actualsize, msg, complainer );
    if ~ok
        if expectedsize > actualsize
            v( (actualsize+1):expectedsize, 1 ) = defaultval;
        elseif expectedsize < actualsize
            v( (expectedsize+1):end ) = [];
        end
    end
end

function [ok,v] = rectifysize2( expectedsize, v, defaultval, msg, complainer )
    actualsize = size(v,1);
    ok = checksize( expectedsize, actualsize, msg, complainer );
    if ~ok
        if expectedsize > actualsize
            v( (actualsize+1):expectedsize, : ) = repmat( defaultval, expectedsize-actualsize, 1 );
        elseif expectedsize < actualsize
            v( (expectedsize+1):end, : ) = [];
        end
    end
end

function ok = checksize( expectedsize, actualsize, msg, complainer )
    ok = expectedsize == actualsize;
    if ~ok
        complainer( ['validmesh:bad', msg], ...
            '%s: Expected %d, found %d.', ...
            msg, expectedsize, actualsize );
    end
end

function result = checkheight( m, field, height, complainer )
    if isfield( m, field ) && (size(m.(field),1) ~= height)
        result = 0;
        complainer( ['validmesh:' field], ...
            'Wrong size of %s: %d, expected %d.', ...
            field, size(m.(field),1), height );
    else
        result = 1;
    end
end


function result = checknumel( m, field, num, complainer )
    if isfield( m, field ) && (numel(m.(field)) ~= num)
        result = 0;
        complainer( ['validmesh:' field], ...
            'Wrong size of %s: %d, expected %d.', ...
            field, size(m.(field),1), num );
    else
        result = 1;
    end
end


function result = isedgeof(m,ei,fi)
%ISEDGEOF(P,E,F,EI,FI)  Test whether edge EI is an edge of cell FI.
    result = 1;
    ep1 = m.edgeends(ei,1);
    ep2 = m.edgeends(ei,2);
    fp1 = m.tricellvxs(fi,1);
    fp2 = m.tricellvxs(fi,2);
    fp3 = m.tricellvxs(fi,3);
    fe1 = m.celledges(fi,1);
    fe2 = m.celledges(fi,2);
    fe3 = m.celledges(fi,3);
    if (fe1==ei)
        if ((ep1~=fp2) || (ep2~=fp3)) && ((ep1~=fp3) || (ep2~=fp2))
            fprintf( 1, 'isedgeof( %d, %d ) fails at fe1\n', ei, fi );
            fprintf( 1, 'ep %d %d fp %d %d %d fe %d %d %d\n', ...
                ep1, ep2, fp1, fp2, fp3, fe1, fe2, fe3 );
            result = 0;
        end
    elseif (fe2==ei)
        if ((ep1~=fp3) || (ep2~=fp1)) && ((ep1~=fp1) || (ep2~=fp3))
            fprintf( 1, 'isedgeof( %d, %d ) fails at fe2\n', ei, fi );
            fprintf( 1, 'ep %d %d fp %d %d %d fe %d %d %d\n', ...
                ep1, ep2, fp1, fp2, fp3, fe1, fe2, fe3 );
            result = 0;
        end
    elseif (fe3==ei)
        if ((ep1~=fp1) || (ep2~=fp2)) && ((ep1~=fp2) || (ep2~=fp1))
            fprintf( 1, 'isedgeof( %d, %d ) fails at fe3\n', ei, fi );
            fprintf( 1, 'ep %d %d fp %d %d %d fe %d %d %d\n', ...
                ep1, ep2, fp1, fp2, fp3, fe1, fe2, fe3 );
            result = 0;
        end
    else
        result = 0;
    end
end

function donothing(varargin)
end
