Let, D=|(1,a,a^2,p),(1,b,b^2,q),(1,c,c^2,r),(1,d,d^2,s)|,
where, p=a^3+bcd, q=b^2+cda, r=c^3+dab, s=d^3+abc.
Note that, D can be split as the sum of two determinants
D_1 and D_2, where,
D_1=|(1,a,a^2,a^3),(1,b,b^2,b^3),(1,c,c^2,c^3),(1,d,d^2,d^3)|, &,
D_2=|(1,a,a^2,bcd),(1,b,b^2,cda),(1,c,c^2,dab),(1,d,d^2,abc)|.
So, D=D_1+D_2.
By R_1xxa, R_2xxb, R_3xxc, R_4xxd, we have,
(abcd)D_2=|(a,a^2,a^3,abcd),(b,b^2,b^3,bcda),(c,c^2,c^3,cdab),(d,d^2,d^3,dabc)|, or,
D_2=|(a,a^2,a^3,1),(b,b^2,b^3,1),(c,c^2,c^3,1),(d,d^2,d^3,1)|.
Next, the interchange between C_1" & "C_4 will give -D_2,
i.e., -D_2=|(1,a^2,a^3,a),(1,b^2,b^3,b),(1,c^2,c^3,c),(1,d^2,d^3,d)|.
Similarly, C_2harr C_4 will yield -(-D_2)=D_2.
:.D_2=|(1,a,a^3,a^2),(1,b,b^3,b^2),(1,c,c^3,c^2),(1,d,d^3,d^2)|.
Finally, C_3harrC_4 will result in -D_2.
:.-D_2=|(1,a,a^2,a^3),(1,b,b^2,b^3),(1,c,c^2,c^3),(1,d,d^2,d^3)|.
But then, this means that, -D_2=D_1, or,
D_1+D_2=D=0, as Respected Maganbhai P. has readily
shown!
color(blue)("Enjoy Maths.!")