When the currents are increasing in the circuit, we have e.m.f.s E_1=-L_1\\dfrac{dI_1}{dt} \\pm M \\dfrac{dI_2}{dt},E_2=-L_2\\dfrac{dI_2}{dt} \\pm M\\dfrac{dI_1}{dt} where \\pm sign appears consistently in both equations and depends on the geometry of the coils and the presence of current.Work done in pushing charges dq_1 and dq_2 through each circuit, respectively, is dW=-E_1dq_1-E_2 dq_2 =L_1\\dfrac{dI_1}{dt}dq_1 \\mp M \\dfrac{dI_2}{dt}dq_1+L_2\\dfrac{dI_2}{dt}dq_2 \\mp M\\dfrac{dI_1}{dt}dq_2 I_1=\\dfrac{dq_1}{dt},dq_1=I_1 dt, I_2=\\dfrac{dq_2}{dt},dq_2=I_2 dt, dW=L_1I_1dI_1+L_2I_2dI_2\\mp (MI_1dI_2+MI_2dI_1) =L_1I_1dI_1+L_2I_2dI_2 \\mp Md(I_1I_2) On integrating the above expression from 0 to final current i.e I_2 ,we have U=\\int dW=L_1 \\int^{I_1}_0I_1dI_1+L_2\\int^{I_2}_0 dI_2 \\mp M \\int^{I_1I_2}_0 d(I_1I_2) =\\dfrac{1}{2}L_1I^2_1 +\\dfrac{1}{2}L_2I^2_2 \\mp MI_1I_2