anchor.diff:          end do
       end function
   
  +    subroutine place(c1, p, id, c2, ok)
  +        integer, intent(in)  :: c1(3, 3, 3), p(:, :), id
  +        integer, intent(out) :: c2(3, 3, 3)
  +        logical, intent(out) :: ok
  +        integer              :: i, x, y, z
  +
  +        c2 = c1
  +        do i = 1, size(p, 2)
  +            x = p(1, i)
  +            y = p(2, i)
  +            z = p(3, i)
  +            if (c1(x, y, z) /= 0) then
  +                ok = .false.
  +                return
  +            end if
  +            c2(x, y, z) = id
  +        end do
  +        ok = .true.
  +    end subroutine
  +
  +    function all_puts(cube, piece, id) result(cubes)
  +        integer, intent(in)                   :: cube(3, 3, 3), piece(:, :), id
  +        type(list)                            :: cubes, rots
  +        integer                               :: i, x, y, z, rmax(3)
  +        integer, dimension(3, size(piece, 2)) :: rot, px, py, pz
  +        integer                               :: c1(3, 3, 3)
  +        logical                               :: ok
  +
  +        cubes%length = 0
  +        !allocate(cubes%s(size(piece, 1), size(piece, 2), 3**3))
  +        allocate(cubes%d(3, 3, 3, 3**3))
  +
  +        rots = all_rots(piece)
  +        do i = 1, rots%length
  +            rot = rots%s(:, :, i)
  +            rmax = maxval(rot, 2)
  +            do x = 0, 3 - rmax(1)
  +                px = rot
  +                px(1, :) = px(1, :) + x
  +                do y = 0, 3 - rmax(2)
  +                    py = px
  +                    py(2, :) = py(2, :) + y
  +                    do z = 0, 3 - rmax(3)
  +                        pz = py
  +                        pz(3, :) = pz(3, :) + z
  +                        call place(cube, pz, id, c1, ok)
  +                        if (ok) then
  +                            cubes%length = cubes%length + 1
  +                            !cubes%s(:, :, cubes%length) = pz
  +                            cubes%d(:, :, :, cubes%length) = c1
  +                        end if
  +                    end do
  +                end do
  +            end do
  +        end do
  +    end function
  +
  +    recursive subroutine search(ps, n, cube, sols)
  +        integer, intent(in)       :: ps(:, :, :), n, cube(3, 3, 3)
  +        type(list), intent(inout) :: sols
  +        type(list)                :: puts
  +        integer                   :: i
  +
  +        if (n == 0) then
  +            sols%length = sols%length + 1
  +            sols%d(:, :, :, sols%length) = cube
  +            return
  +        end if
  +
  +        ! TODO Remove searchs with tiny islands.
  +        puts = all_puts(cube, ps(:, :, n), n)
  +        ! TODO Remove puts that are rotations of each other.
  +        do i = 1, puts%length
  +            call search(ps, n - 1, puts%d(:, :, :, i), sols)
  +        end do
  +    end subroutine
  +
       subroutine read_pieces(vecs)
           integer, allocatable, intent(inout) :: vecs(:, :, :)
           character(len=20) :: name_
  
positive.code:                 do k = 0, 3
                      p1 = push_to_one(p1)
                      d = sparse_to_dense(p1)
                      if (.not.contains_p(ps%d, p1)) then
                          ps%length = ps%length + 1
                          ps%val(:, :, ps%length) = p1
                          ps%d(:, :, :, ps%length) = d
                      end if
                      p1 = rotate(p1, 1, 1)
  
negative.code: 
      integer, parameter   :: isin(0:3) = (/ 0, 1,  0, -1 /), &
                              icos(0:3) = (/ 1, 0, -1,  0 /)
      integer              :: mrot(3, 3, 3, 4), cube(3, 3, 3)
      integer, allocatable :: vecs(:, :, :)
      type(list)           :: rots, sols, test
      integer              :: ii
  
      cube = 0
      call init_mrot(mrot)
      call read_pieces(vecs)
      call print_pieces(vecs)
      sols%length = 0
      allocate(sols%d(3, 3, 3, 128))
  
      call search(vecs, size(vecs, 3), cube, sols)
      print *, sols%length
      do ii = 1, sols%length
          call print_pieces(sols%d(:, :, :, ii))
          print *, "********************************"
      end do
  
  contains
      subroutine init_mrot(mrot)
          integer, intent(inout) :: mrot(3, 3, 3, 4)
  