The induced electric and magnetic surface currents on each point of the scatterer object can be calculated from the Magnetic and Electric Field Integral Equations (MFIE & EFIE):
:<math> \mathbf{J(r)} = (1+\alpha)\mathbf{\hat{n}} \times \left\lbrace \begin{align} & \mathbf{ H^{inc}(r) } - jk_0 \iint_{S_J} \left( 1 - \frac{j}{k_0 R} \right) (\mathbf{ \hat{R} \times J(r') }) \frac{e^{-jk_0 R}}{4\pi R} \,ds' \\ & -j k_0 Y_0 \iint_{S_M} \left[ \left( 1 - \frac{j}{k_0 R} - \frac{1}{(k_0 R)^2} \right) \mathbf{M(r')} - \left( 1-\frac{3j}{k_0 R}-\frac{3}{(k_0 R)^2} \right) \mathbf{ (\hat{R} \cdot M(r')) \hat{R} } \right] \frac{e^{-jk_0 R}}{4\pi R} \,ds' \end{align} \right\rbrace </math>
:<math> \mathbf{M(r)} = -(1-\alpha)\mathbf{\hat{n}} \times \left\lbrace \begin{align} & \mathbf{ E^{inc}(r) } + jk_0 \iint_{S_M} \left( 1 - \frac{j}{k_0 R} \right) (\mathbf{ \hat{R} \times M(r') }) \frac{e^{-jk_0 R}}{4\pi R} \,ds' \\ & -j k_0 Z_0 \iint_{S_J} \left[ \left( 1 - \frac{j}{k_0 R} - \frac{1}{(k_0 R)^2} \right) \mathbf{J(r')} - \left( 1-\frac{3j}{k_0 R}-\frac{3}{(k_0 R)^2} \right) \mathbf{ (\hat{R} \cdot J(r')) \hat{R} } \right] \frac{e^{-jk_0 R}}{4\pi R} \,ds' \end{align} \right\rbrace </math>
:<math> \mathbf{J(r)} = \mathbf{\hat{n}\times} \begin{bmatrix} 1-R_{||} & 0 \\ 0 & 1-R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{H_{||}^{tot}(r)} \\ \mathbf{H_{\perp}^{tot}(r)} \end{bmatrix} </math>
:<math> \mathbf{M(r)} = -\mathbf{\hat{n}\times} \begin{bmatrix} 1+R_{||} & 0 \\ 0 & 1+R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{E_{||}^{tot}(r)} \\ \mathbf{E_{\perp}^{tot}(r)} \end{bmatrix} </math>
<!--[[File:PO9(1).png]]--> where '''R''' ='''r''' - '''r'''', R = |'''R'''|, and :<math>\mathbf{ \hat{R} = \frac{R}{|R|} = \frac{r-r'}{|r-r'|} }</math><!--[[File:PO11.png]]--> The shadowing phenomenon can indeed be attributed to near-field interaction of surface currents. The current on the lit region produces a scattered field in the forward direction that is almost equal and out of phase with the incident wave. Hence, the sum of the scattered field and incident field over the shadowed region almost cancel each other, giving rise to a very small field there. This suggests that keeping track of multiple scattering can take care of shadowing problems automatically. In addition, the effects of multiple scattering can be readily accounted for by an iterative PO approach to be formulated next.
The starting point for the iterative PO solution is the above MFIE and EFIE integral equations. To the first (zero-order) approximation, we can write
:<math> \begin{align} & \mathbf{J^{(0)}(r)} = (1+\alpha) \mathbf{ \hat{n} \times H^{inc}(r) } \\ & \mathbf{M^{(0)}(r)} = -(1-\alpha) \mathbf{ \hat{n} \times E^{inc}(r) } \end{align} </math>
<!--[[File:PO13.png]]-->
:<math> \mathbf{J^{(0)}(r)} = \mathbf{\hat{n}\times} \begin{bmatrix} 1-R_{||} & 0 \\ 0 & 1-R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{H_{||}^{inc}(r)} \\ \mathbf{H_{\perp}^{inc}(r)} \end{bmatrix} </math>
:<math> \mathbf{M^{(0)}(r)} = -\mathbf{\hat{n}\times} \begin{bmatrix} 1+R_{||} & 0 \\ 0 & 1+R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{E_{||}^{inc}(r)} \\ \mathbf{E_{\perp}^{inc}(r)} \end{bmatrix} </math>
which are the conventional PO currents. However, this approximation does not formally recognize the lit and shadowed areas. Instead of identifying the exact boundaries of the lit and shadowed areas over a complex target, a simple condition is used first to find the primary shadowed areas. Then, through PO iterations all shadowed areas are determined automatically. When calculating the field on the scatterer for every source point, a primary shadowing condition given by '''n.k'''< 0 is examined. In complex scatterer geometries, there are shadowed points in concave regions where '''n.k'''> 0, but the correct shadowing is eventually achieved through the iteration of the PO currents. Therefore, in computation of the above equations, only the contribution of the points that satisfy the following condition are considered:
At the subsequent iterations, the higher order PO currents are given by:
:<math> \mathbf{J^{(n)}(r)} = (1+\alpha)\mathbf{\hat{n}} \times \left\lbrace \begin{align} & - jk_0 \iint_{S_J} \left( 1 - \frac{j}{k_0 R} \right) (\mathbf{ \hat{R} \times J^{(n-1)}(r') }) \frac{e^{-jk_0 R}}{4\pi R} \,ds' \\ & -j k_0 Y_0 \iint_{S_M} \left[ \left( 1 - \frac{j}{k_0 R} - \frac{1}{(k_0 R)^2} \right) \mathbf{M^{(n-1)}(r')} - \left( 1-\frac{3j}{k_0 R}-\frac{3}{(k_0 R)^2} \right) \mathbf{ (\hat{R} \cdot M^{(n-1)}(r')) \hat{R} } \right] \frac{e^{-jk_0 R}}{4\pi R} \,ds' \end{align} \right\rbrace </math>
:<math> \mathbf{M^{(n)}(r)} = -(1-\alpha)\mathbf{\hat{n}} \times \left\lbrace \begin{align} & jk_0 \iint_{S_M} \left( 1 - \frac{j}{k_0 R} \right) (\mathbf{ \hat{R} \times M^{(n-1)}(r') }) \frac{e^{-jk_0 R}}{4\pi R} \,ds' \\ & -j k_0 Z_0 \iint_{S_J} \left[ \left( 1 - \frac{j}{k_0 R} - \frac{1}{(k_0 R)^2} \right) \mathbf{J^{(n-1)}(r')} - \left( 1-\frac{3j}{k_0 R}-\frac{3}{(k_0 R)^2} \right) \mathbf{ (\hat{R} \cdot J^{(n-1)}(r')) \hat{R} } \right] \frac{e^{-jk_0 R}}{4\pi R} \,ds' \end{align} \right\rbrace </math>
:<math> \mathbf{J^{(n)}(r)} = \mathbf{\hat{n}\times} \begin{bmatrix} 1-R_{||} & 0 \\ 0 & 1-R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{H_{||}^{(n-1)}(r)} \\ \mathbf{H_{\perp}^{(n-1)}(r)} \end{bmatrix} </math>
:<math> \mathbf{M^{(n)}(r)} = -\mathbf{\hat{n}\times} \begin{bmatrix} 1+R_{||} & 0 \\ 0 & 1+R_{\perp} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{E_{||}^{(n-1)}(r)} \\ \mathbf{E_{\perp}^{(n-1)}(r)} \end{bmatrix} </math>
<!--[[File:PO9(1).png]]-->
<!--[[File:PO14(1).png]]-->
For most practical applications, iterations up to the second order is sufficient. The iterative solution will not only account for double-bounce scattering over the lit regions but it also removes the lower order currents erroneously placed over concave shadowed areas.
