The Cartesian parallel magnetic imaging problem is formulated variationally using a high-order penalty for coil sensitivities and a total variation like penalty for the reconstructed image. Then the optimality system is derived and numerically discretized. The objective function used is non-convex, but it possesses a bilinear structure that allows the ambiguity among solutions to be resolved technically by regularization and practically by normalizing a pre-estimated norm of the reconstructed image. Since the objective function is convex in each single argument, convex analysis is used to formulate the optimality condition for the image in terms of a primaldual system. To solve the optimality system, a nonlinear GaussSeidel outer iteration is used in which the objective function is minimized with respect to one variable after the other using an inner generalized Newton iteration. Computational results for in vivo MR imaging data show that a significant improvement in reconstruction quality can be obtained by using the proposed regularization methods in relation to alternative approaches.