(*Remark: Not all the presented options and parameters used for the calculation are necessary. There are probably more efficient methods to implement such derivation in WM14. Therefore, we present here what worked for us, without aiming to optimize the derivations.*)
Clear["Global`*"];

(*allow more assumptions*)
SetSystemOptions["SimplificationOptions"->{"AssumptionsMaxNonlinearVariables"->25,"AssumptionsMaxVariables"->25}];

(*Assumptions*)
$Assumptions=bx>0 && by>0&&hx>0 && hy>0&& Lx>0&& q2>0&&  q1 ∈ Reals && bx+q2+2*by>=y1>=bx+q2&&bx+q2>x1 &&Element[m,Integers]&&m>3&&v>1&&p>0 && Element[x1,Reals]&& Element[x2,Reals]&& Element[y1,Reals]&& Element[y2,Reals]&&y1>x1;

(*A function that helps integration*)
integrateD[x_,v_,l_,u_]:=Module[{int=x,var=v,low=l,up=u},
indef=FullSimplify[Integrate[int,var],TimeConstraint->1];
upper=Simplify[Limit[indef,var->up,Direction->"FromBelow"],TimeConstraint->1];
lower=Simplify[Limit[indef,var->low,Direction->"FromAbove"],TimeConstraint->1];
res=FullSimplify[upper-lower,TimeConstraint->1];
Return[res]]

(*We start from a point-pair interaction of (m^th) order, 1/(r^m). The coordinates of the first particle are (x1,x2,0) and of the second (y1,y2,q1). Therefore, the distance between them is r=(Sqrt[(p^2)+(q1^2)]), where p is the in-plane distance. For Cartesian coordinates p=(Sqrt[(((y1-x1))^2)+(((y2-x2))^2)]). If we integrate along q1, we get a point-infinite line (along the e3-axis) law, valid for m>1.*)
pointInfLine3=FullSimplify[Integrate[(1/(((p^2+(q1^2)))^(m/2))),{q1,-Infinity,Infinity}]]

(*This integration is the first step when considering the interaction between a finite-sized body and an infinite body, see e.g. Langbein 1972 for m=6. This result suggests that the (m^th) order interaction potetnial between a point and an infinite line along the e3 axis can be represented as a factor fm=(((Sqrt[π]) Gamma[(1/2) (-1+m)])/(Gamma[(m/2)])) multiplied by the point-pair interaction potential of (((m-1))^th) order, (1/(p^(m-1))). We can generalize this: an interaction of (m^th) order between in-plane points/figures can be used to find an interaction between the same objects, of which one is infinitely extended perpendicular to the plane - RIIPI (the rule of infinite and in-plane body interaction), cf. paper.*)

(*Remark: let us note that a point-finite line interaction law is much more complicated and involves a hypergeometric function. For odd m, it evaluates to rational functions, while for even m it consists of a hypergeometric function*)
pointFinLine3=FullSimplify[Integrate[(1/(((p^2+(q1^2)))^(m/2))),{q1,-Lx,Lx}]]
FullSimplify[pointFinLine3/.m->{3,4,5,6,7,8}]

(*Next, let us derive a point-infinite line (along e2 axis) in-plane law. The point has coordinates (x1,x2) while the line is on coordinate y1. Valid for m>1*)
pointInfLine2=Integrate[1/((y1-x1)^2+(y2-x2)^2)^(m/2),{y2,-Infinity,Infinity}]

(*The result is analogous to the integration along the e3 axis, as expected. Therefore, now we can get a point-infinite plane law either from the P-L2 law or from the P-L3 law, by simply using RIIPI*)
fm=(((Sqrt[π]) Gamma[(1/2) (-1+m)])/(Gamma[(m/2)]));
pointInfPlane23=FullSimplify[fm*(pointInfLine2/.{m->m-1})]
(*or*)
pointInfPlane23=FullSimplify[fm*(pointInfLine3/.{m->m-1,p->y1-x1})]

(*Next, let us find a point-infinite plate (thickness 2by) and a point-infinite half-space laws, by integrating the point-infinite plane law along the coordinate y1*)
pointInfPlate23=FullSimplify[Integrate[pointInfPlane23,{y1,x1+q2,x1+q2+2by}]]
(*manual simplification*)
pointInfPlate23=((2 π (-(((2 by+q2))^(3-m))+(q2^(3-m)) ))/((-3+m) (-2+m)));
FullSimplify[((2 π (((q2 (2 by+q2)))^(-m)) (-(q2^m) (((2 by+q2))^3)+(q2^3) (((2 by+q2))^m)))/((-3+m) (-2+m)))-pointInfPlate23]
(*The point-HS law can be obtained either by integrating point-infinite plane from x1+q2 to Infinity, or by letting by->Infinity in the point-PT law. m>4*)
pointInfHS23=FullSimplify[Integrate[pointInfPlane23,{y1,x1+q2,Infinity}]]
FullSimplify[Limit[pointInfPlate23,by->Infinity,Assumptions->m>3]]

(*The interaction potential point-infinite plane (parallel to e2e3 plane) is also an interaction potential between a unit-area plane figure and an infinite plane. By multiplying it with 2Lx*2hx, we obtain a potential between a plane figure (sheet/rectangle) with dimensions 2Lx/2hx and an infinite plane*)
planeInfPlane23=2hx 2 Lx pointInfPlane23

(*By integrating point-infinite plane law over x1, we obtain a unit-area plate with thickness 2bx (or a line of length 2bx along the e1 axis)-infinite plane law, where y1 is the distance between the middle of the plate/line and the infinite plane. Valid for m>3*)
uPlateInfPlane23=FullSimplify[Integrate[pointInfPlane23,{x1,-bx,bx}]]

(*manual simplification*)
uPlateInfPlane23=((2 π (-(((bx+y1))^(3-m))+(((-bx+y1))^(3-m))))/((-3+m) (-2+m)))
FullSimplify[uPlateInfPlane23-((2 π (((-(bx^2)+(y1^2)))^(-m)) (-(((-bx+y1))^m) (((bx+y1))^3)-(((bx-y1))^3) (((bx+y1))^m)))/((-3+m) (-2+m)))]

(*Finally, by integrating over the thickness, along y1, we get a plate (thickness 2bx, unit area)-infinite plate (thickness 2 by) law, where q2 is the gap between plates. Valid for m>4*)
plateInfPlate23=FullSimplify[integrateD[uPlateInfPlane23,y1,q2+bx,q2+bx+2 by]]

(*vdW case, m=6*)
plateInfPlate23/.m->6

(*If we let by->Infinity, it is a unit-area plate-infinite HS law*)
plateHS=Limit[plateInfPlate23,by->Infinity,Assumptions->m>4]

(*Finally, by letting bx->Infinity, it is a unit-area HS-infinite HS law. We can find it in Israelachvili p. 211*)
hsHS[q2_]=Limit[plateHS,bx->Infinity,Assumptions->m>4]

(*The rule of infinite half-body interaction (RIHBI) is evident from the plate-plate and HS-HS interactions. In order to find PT-PT from HS-HS, we need to add HS-HS(q2) and HS-HS(q2+2bx+2by) and to subtract HS-HS(2bx+q2) and HS-HS (2by+q2)*)

FullSimplify[hsHS[q2]+hsHS[q2+2bx+2by]-hsHS[q2+2bx]-hsHS[q2+2by]]

(*We can also follow another order of integration to get slightly different laws. We are especially interested in section-section, section-plate, and section-HS laws.
An in-plane rectangle-rectangle law is derived in Notebook 5. Here, we present rectangle-infinite plane 23, rectangle-infinite plate 23, and rectangle-infinite HS laws. The rectangle is in the e2e1 plane and has dimensions 2bx/2hx*)
rectSecInfPlane23=uPlateInfPlane23*2*hx
rectSecInfPlate23=plateInfPlate23*2*hx
rectSecHS23=plateHS*2*hx
(*From these expressions, we can obtain rectangular prism-infinite plane 23, rectangular prism-infinite plate 23, and rectangular prism-infinite HS 23 laws. The rectangular prism has a length of 2*Lx (along the e3 direction), while the cross-sectional dimensions are 2bx/2hx)*)
rectPrismInfPlane23=rectSecInfPlane23*2*Lx
rectPrismInfPlate23=rectSecInfPlate23*2*Lx
rectPrismHS23=rectSecHS23*2*Lx