The Dirac electron in graphene placed in magnetic fields orthogonal to the layer is studied. The determination of the eigenvalues makes use of two auxiliary Schrödinger Hamiltonians intertwined by a first-order differential operator. The magnetic field is initially chosen such that the two auxiliary potentials are shape invariant. Then, more general Darboux transformations are used, producing new magnetic fields for which there are again analytic solutions to the graphene problem. The iterations of the method in the real, singular, and complex case are discussed, trying to produce Hermitian graphene Hamiltonians with applied magnetic fields. The case with periodic and quasiperiodic magnetic fields is as well analyzed.