Two approximations of the Hessian matrix as limited-memory operators are built from the limited-memory BFGS inverse Hessian approximation provided by the minimization code, in view of the specification of the inverse analysis/forecast error covariance matrix in variational data assimilation. Some numerical experiments and theoretical considerations lead to reject the limited-memory DFP Hessian approximation and to retain the BFGS one for the applications foreseen. Conditioning issues are explored and a preconditioning strategy via a change of control variable is proposed, based on a suitable Cholesky factorization of the limited-memory inverse Hessian matrix. This factorization is implemented as the composition of linear operators. The memory requirements and the number of floating-point operations required by the method are given and confirmed by numerical experiments. The method is found to have a strong potential for variational data assimilation systems using high resolution ocean or atmosphere general circulation models.