As promised, Lagrange published his own method in the next issue of Miscellanea Taurinensia [11, p. 179]. His “De différens problèmes de calcul intégral” (“On different problems on integral calculus”) immediately preceded d'Alembert's paper in the journal (see Figure 7). In order to make his method more understandable, we'll adopt the notation proposed by historian S. S. Demidov in [4, p. 370].

On the integration of the equation
\[Ly + M \frac{dy}{dt}+N\frac{d^2y}{dt^2}+P\frac{d^3y}{dt^3}+\&c. = T \hspace{1 cm}(A)\]
in which \(L,M,N \text{ etc. } T\) are functions of \(t\).
1. I multiply this equation by \(\xi dt\), \(\xi\) being an undetermined variable, I take the integral, I have
\[\int L \xi y dt + \int M \xi \frac{dy}{dt} dt+ \int N \xi \frac{d^2y}{dt^2} dt + \int P \xi \frac{d^3y}{dt^3}dt + \&c. = \int T \xi dt ;\]
I change the expressions \(\int M \xi \frac{dy}{dt} dt, \quad \int N \xi \frac{d^2y}{dt^2} dt ,\quad \int P \xi \frac{d^3y}{dt^3} dt\) to their equivalents \(M \xi y  \int \frac{d M \xi}{dt} y\,dt,\) \(N \xi \frac{dy}{dt}  \frac{d N \xi}{dt} y + \int \frac{d^2 N \xi}{dt^2} y\,dt,\) \(P \xi \frac{d^2y}{dt^2}\frac{d P \xi}{dt} \frac{dy}{dt}+ \frac{d^2 P \xi}{dt^2} y  \int \frac{d^3 P \xi}{dt^3} y\,dt\,\,\&c.\) And I have, by ordering the terms by relation to \(y,\) ....

Figure 7. Lagrange began his article on his method for reduction of order [11, p. 179]. Lagrange's French is followed by the authors' translation into English.
Suppose we are given a linear differential equation of the form
\[L_m(y) = a_m(x) y^{(m)} +a_{m1}(x) y^{(m1)} + \dots + a_0(x) y = X(x)\]
and we know that \(y_1\) is a solution to the associated homogenous equation \(L_m(y)=0\).
Lagrange multiplied both sides of this differential equation by an undefined function \(z\) and integrated, giving
\[\int z\, [a_m(x) y^{(m)} +a_{m1}(x) y^{(m1)} + \dots + a_0(x) y ] dx = \int z X(x) dx.\] He executed \(k\) integrations by parts on each \(z a_{k}(x) y^{(k)},\) which left \((z a_{k}(x))^{(k)} y\) under the integral. For example when \(k=2\), \[\int z a_{2}(x) y'' dx = z a_{2}(x) y' \int (z a_{2}(x))' y' dx\]
\[= z a_{2}(x) y'(z a_{2}(x))' y+\int (z a_{2}(x))'' y \, dx.\] Doing this for each term on the left side and collecting terms under the integral gave
\[\int z L_n(y) dx = A(x,y,y', \dots, y^{m1}, z, z', \dots, z^{(m1)}) + \int y L_m^*(z) dx\] \[= \int z X(x) dx.\] The expression \(A(x,y,y', \dots, y^{m1}, z, z', \dots, z^{(m1)})\) is called the bilinear concomitant and is a bilinear function in \(y\) and \(z\) and their derivatives. (In general a concomitant is “a phenomenon that naturally accompanies or follows something.”) Finally, the adjoint \(L_m^*(z)\) of the differential operator \(L_n(y)\) is given by
\[L_m^*(z) =(1)^m (za_m(x))^{(m)}+(1)^{(m1)}(z a_{m1}(x))^{(m1)} +\dots  (z a_{1}(x) )' + z a_0(x).\]
If \(L_m=L_m^*\), then the operator is selfadjoint [16, p. 96], in which case any solution \(y_1=z_1\) of \(L_m(y)\) is also a solution to \(L_m^*(z)\). It follows that the integral \(\int y L_m^*(z_1)dx =0\) and the equation
\[A(x,y,y', \dots, y^{m1}, z, z', \dots, z^{(m1)}) + \int y L_m^*(z) dx = \int z X(x) dx\] would simplify to the \((m1)\)th order linear equation
\[A(x,y,y', \dots, y^{m1}, z_1, z_1', \dots, z_1^{(m1)}) + 0 = \int z_1 X(x) dx,\] as desired.
Further reduction can be accomplished when additional solutions are known. If \(L\) is selfadjoint and \(z_1, z_2, \dots, z_k\) are \(k\) solutions to \(L_m^*(z)=0\) (since \(k\) solutions \(y_1, y_2, \dots, y_k\) are known to \(L_m(y)=0\)), we obtain a system of \(k\) linear differential equations
\(A(x,y,y', \dots, y^{m1}, z_1, z_1', \dots, z_1^{(m1)}) + 0\) 
\(= \int z_1 X(x)\) 
\(A(x,y,y', \dots, y^{m1}, z_2, z_2', \dots, z_2^{(m1)}) + 0\) 
\(= \int z_2 X(x)\) 

\(\vdots\) 
\(A(x,y,y', \dots, y^{m1}, z_k, z_k', \dots, z_k^{(m1)}) + 0\) 
\(= \int z_k X(x)\) 
and \(y^{(m1)}, y^{(m2)}, \dots, y^{(mk+1)}\) can be eliminated to obtain a differential equation of order \((mk)\) [15, p. 124].
Certainly not every operator is selfadjoint, though those that describe physical processes often are. Finding conditions for an operator to be selfadjoint is easy. By equating the coefficients in the differential equation \[L_m(y) = a_m(x) y^{(m)} +a_{m1}(x) y^{(m1)} + \dots + a_0(x) y = X(x)\] and with those of the adjoint \[L_m^*(z) =(1)^m (za_m(x))^{(m)}+(1)^{(m1)}(z a_{m1}(x))^{(m1)} +\dots  (z a_{1}(x) )' + z a_0(x),\] the following theorem arises:
Theorem [16, p. 98]. The secondorder linear differential equation
\[a_2(x)y''+a_1(x)y'+a_0(x)y = 0\] is selfadjoint if \(a_1(x)=a_2'(x)\) and the thirdorder linear differential equation
\[a_3(x)y'''+a_2(x)y''+a_1(x)y'+a_0(x)y=0\]
is selfadjoint if \(a_2(x)=\frac{3}{2}a_3'(x)\) and \(a_0(x)=\frac{1}{2}(a_1(x)\frac{1}{3}a_2'(x))'\).
When the differential equation is not selfadjoint, Lagrange repeated the process on \(L_m^{*}(z)=0\), or equivalently
\[(1)^m (za_m(x))^{(m)}+(1)^{(m1)}(z a_{m1}(x))^{(m1)} +\dots  (z a_{1}(x) )' +z a_0(x) =0.\] He multiplied both sides by an unknown function \(w(x)\) and integrated. Executing multiple integrations by parts gives
\[\int w L_m^{*}(z) dx = B(x,z,z',\dots,z^{(m1)}, w, w', \dots, w^{(m1)}) + \int z L_m^{**}(w) dx\]
\[= \int z \cdot 0 \, dx = c,\] where \(L_m^{**}\) is the adjoint of the adjoint of the original operator. Lagrange proved that \(L_m^{**}=L_m\). Hence a solution, \(w_1\), to the differential equation \(L_m^{**}(w) =0\) is known; namely \(w_1=y_1\). Plugging this into the preceding equations shows that
\[B(x,z,z',\dots,z^{(m1)}, w_1, w_1', \dots, w_1^{(m1)}) =c,\] which is a differential equation of order \(m1\). If we integrate this equation to find a solution, \(z_1\), then we have that \(L_m^{*}(z_1)=0\). Finally, by the equation \[A(x,y,y', \dots, y^{m1}, z, z', \dots, z^{(m1)}) + \int y L_m^*(z) dx = \int z X(x) dx,\] we have \[A(x,y, y', \dots, y^{(m1)}, z_1,z_1',\dots,z_1^{(m1)}) = \int z_1 X(x) dx,\]
which is also a differential equation of order \(m1\).