AdaptFromFlags_QuadTreeMesh2D Subroutine

  subroutine AdaptFromFlags_QuadTreeMesh2D(this,flag)
    !! Mutate the forest from a per-leaf flag array (indexed over the current leaves in this%leaf
    !! order, e.g. the flag array produced by SELF_RefinementIndicator_2D):
    !!
    !!   flag(i) = QUADTREE_REFINE  (+1) -> subdivide leaf i
    !!   flag(i) = QUADTREE_COARSEN (-1) -> merge leaf i's family if ALL four siblings are leaves
    !!                                      and ALL are flagged COARSEN (standard coarsening rule)
    !!   flag(i) = QUADTREE_KEEP     (0) -> unchanged
    !!
    !! Refinement and coarsening are both resolved against the pre-adaptation leaf snapshot, so
    !! the two operations never interfere. The leaf list is rebuilt on return.
    implicit none
    class(QuadTreeMesh2D),intent(inout) :: this
    integer,intent(in) :: flag(1:this%nLeaves)
    ! Local
    integer :: i,c,p,node,nSnap
    integer,allocatable :: flagOfNode(:)
    integer,allocatable :: refineList(:)
    logical :: family

    nSnap = this%nLeaves

    ! Map the flag onto node ids so coarsening can test whole families by node.
    allocate(flagOfNode(1:this%nNodes))
    flagOfNode = QUADTREE_KEEP
    allocate(refineList(1:nSnap))
    refineList = 0
    do i = 1,nSnap
      flagOfNode(this%leaf(i)) = flag(i)
      if(flag(i) == QUADTREE_REFINE) refineList(i) = this%leaf(i)
    enddo

    ! ---- Coarsening : merge families all of whose four leaf children are flagged COARSEN ----
    do i = 1,nSnap
      if(flag(i) /= QUADTREE_COARSEN) cycle
      node = this%leaf(i)
      p = this%parent(node)
      if(p == 0) cycle ! a root cannot be coarsened
      ! Only act once per family (when processing its first child), and only if every child is a
      ! leaf flagged COARSEN.
      if(this%child(1,p) /= node) cycle
      family = .true.
      do c = 1,4
        if(this%child(1,this%child(c,p)) /= 0) family = .false. ! child not a leaf
        if(flagOfNode(this%child(c,p)) /= QUADTREE_COARSEN) family = .false.
      enddo
      if(family) this%child(1:4,p) = 0 ! detach children -> p becomes a leaf again
    enddo

    ! ---- Refinement : subdivide flagged leaves that are still leaves ----
    do i = 1,nSnap
      if(refineList(i) == 0) cycle
      node = refineList(i)
      if(this%child(1,node) == 0) call this%RefineNode(node)
    enddo

    deallocate(flagOfNode,refineList)

    call this%RebuildLeaves()

  endsubroutine AdaptFromFlags_QuadTreeMesh2D