Implementation

For the driven cavity problem, discretized by 3 families of parallel lines, the boolean function is simple.

int p2node(int k)
{     // returns true if vertex k is a velocity only node, 
    int j = k/(nx+1);
    int i = k - (nx + 1) * j;
    return (2*(int(i/2)) != i)||(2*(int(j/2)) != j);
}

The function which sets the pressure at the P2 nodes to their correct value uses the edges:

void p2pressure(Grid& g, Vector<vectNB>& f)
{// assumes pressure at p2 vertices is zero
    for(int k =0; k< g.ne; k++)
    {    if(p2node(g.no(g.e[k].in ))&& !p2node(g.no(g.e[k].out ))) 
            f[g.no(g.e[k].in )][2] += f[g.no(g.e[k].out)][2]/2;
        if(p2node(g.no(g.e[k].out))&& !p2node(g.no(g.e[k].in ))) 
            f[g.no(g.e[k].out)][2] += f[g.no(g.e[k].in )][2]/2;
    }
}

Finally the function which build the P1-iso-P2/P1 matrix from the P1/P1 matrix loops also on the edges to find the nonzero entries of the P1/P1 matrix and then find the neighbors at their neighbors by the neighbor arrays g.v[i].mate[].

void fromp1top2(Grid& g, Bandmatrix<matNB,vectNB>& aa)
{
  int i,i1,j,j1,k,m,l, ii[2];
  for( i1=0; i1<g.ne; i1++)
  { 
      ii[0] = g.no(g.e[i1].in); ii[1] = g.no(g.e[i1].out);
      for(l=0;l<2;l++)
      {    i = ii[l==1?0:1];
          if(k=ii[l], p2node(k))
            for( j1=0; j1<g.v[k].nmate; j1++)
                if(j = g.no(g.v[k].mate[j1]),!p2node(j)) 
                    for(m=0;m<2;m++) 
                    {
                        aa(i,j)(m,2) += aa(i,k)(m,2)/2;
                        aa(j,i)(2,m) += aa(k,i)(2,m)/2;
                    }
        if(p2node(i)) aa(i,i)(2,2) = 1e10;
      }
  }        
}