Variation of parameters and solutions of composite products of linear differential equations

Given a basis of solutions to k ordinary linear differential equations lj[y]=0 (j=1,2,…,k), we show how Green's functions can be used to construct a basis of solutions to the homogeneous differential equation l[y]=0, where l is the composite product l=l1l2…lk. The construction of these solutions is elementary and classical. In particular, we consider the special case when l=l1k. Remarkably, in this case, if {y1,y2,…,yn} is a basis of l1[y]=0, then our method produces a basis of l1k[y]=0 for any k∈N. We illustrate our results with several classical differential equations and their special function solutions.