Computing derivatives of observables with respect to parameters of the
theory is a powerful tool in lattice QCD, as it allows the study of physical
effects not directly accessible in the original Monte Carlo simulation.
Prominent examples of this include the impact of the up-down quark mass
difference and electromagnetic corrections. In this work, we present a new
approach based on automatic differentiation to evaluate such derivatives to
arbitrarily high orders, where particular emphasis will be placed on strong
isospin-breaking effects and on the propagation of derivatives through the
conjugate gradient algorithm in the computation of correlation functions.