  An algorithm for imaginary time evolution of a fermionic projected entangled
pair state (PEPS) with ancillas from infinite temperature down to a finite
temperature state is presented. As a benchmark application, it is applied to
spinless fermions hopping on a square lattice subject to $p$-wave pairing
interactions. With a tiny bias it allows to evolve the system across a
high-temperature continuous symmetry-breaking phase transition.
