  Higher order calculations in perturbative Quantum Field Theories often
produce coupled linear systems of differential equations which factorize to
first order. Here we present an algorithm to solve such systems in terms of
iterated integrals over an alphabet the structure of which is implied by the
coefficient matrix of the given system. We apply this method to calculate the
master integrals in the color-planar and complete light quark contributions to
the three-loop massive form factors.
