1. Introduction

The angular spectrum method is a rigorous method to synthesize near and far-field electromagnetic (EM) and acoustic fields based on the principles of diffraction theory. Its underlying concept stems from the fact that any of these fields can be presented as a superposition of individual plane waves, each with its unique amplitude, phase, and direction [1–6]. The ASM provides a robust framework for examining wave propagation in various scenarios by decomposing the fields into constituent plane waves.

The ASM has two main advantages in electromagnetic simulation and analysis [7–16]. First, due to its source-free eigenfunction nature, the incident field can be synthesized from known planar electric field distribution without singularities. Secondly, it addresses complex wave phenomena, including diffraction, interference, and scattering, often encountered in practical applications [17]. However, ASM is limited to synthesizing EM fields from planar field distributions. The difficulty with handling curved, non-planar surfaces hampers its applicability in various fields that require precise control and manipulation of electromagnetic beams, such as optics, imaging, and beam-shaping applications.

This article proposes a curved boundary integral method (CBIM), which is closely related to the 2D ASM, except that the surface can be curved in the CBIM. Both methods are locally equal on the propagation axis’s positive side but have some exceptions on the negative side. Thus CBIM is strictly a source method, while 2D ASM is a source-free method, allowing the beams to propagate through a field distribution.

In [18], the planar ASM was expanded to approximate diffractive fields from 2D curved surfaces. The expansion was obtained by dividing the curved 2D surface into the stepwise subregions windowed by a Gaussian function. This approach was later expanded to the 3D surfaces in [17], where the authors used the same planar stepwise subsections with Gaussian distribution. Both methods approximate well the diffraction, reflections, and transmission from the curved single-layer boundaries when the radius of curvature (RoC) is much larger than the wavelength (e.g. >100) [19].

The CBIM derivation begins with the same approach of dividing an arbitrary surface into planar subregions, presented as equal amplitude areas without Gaussian distribution or windowing (e.g. Stepwise AS). Then the areas of these subregions are shrunk into infinitesimally small source points, which approach a Dirac delta function in the limit. This approach yields exact results (in the limit of infinitesimally small partitions or infinite frequency) of the diffractive field resulting from an arbitrary surface distribution considering surface structures comparable to or larger than the illumination wavelength, i.e., $k\alpha \geq 4\pi$, where $k$ is the wavenumber and $\alpha$ is the radius of the sphere. Addtitionally, the advantages of the presented method are non-singularity at the source points, allowing the synthesis from any compactly supported and continuous electric field distribution on a regular surface. We have utilized these results in [20–23].

Moreover, each source point is analogous to the angular spectrum presentation of the Stratton-Chu method with added polarization components. The Stratton-Chu method obtains the electromagnetic radiation from a source point as the curl of a magnetic dipole multiplied by scalar Green’s function [24]. This presentation satisfies Maxwell’s equations and approximates well-curved surfaces when the RoC is larger than a few wavelengths [25].

2. Curved boundary integral method

Section 2.1 proposes a CBIM presentation for any compactly supported and continuous electric field distribution defined on a regular surface. Otherwise the distribution can be chosen arbitrarily. Section begins by introducing a general parameterization of an arbitrary surface. Then the local base vectors are derived from parametrization to determine a Cartesian coordinate system for surface sources. Further, the electric field polarization is modified by base vector rotation. Also, a transformation matrix is presented for mapping the local source coordinates in a global coordinate system and vice versa. After coordinate mapping, an arbitrary surface parameterization is discretized into infinitesimally small surface elements. The key concept is that Riemann’s surface integral combines differential surface elements [26], which can be approximated as locally planar elements. The CBIM models the radiated field from each element, and the total electromagnetic field is the sum of the modeled fields from each source point by the superposition principle [1]. This approach relies on the assumption that when the areas of the elements are sufficiently small, the synthesized field from the differential elements approaches the field from the original surface. Section 2.2 details a mathematical proof of the CBIM presentation. Section 2.3 expands the CBIM presentation for the magnetic field. Section 2.4 introduces the repropagation method between surfaces.

2.1 Electric field representation

Let the global Cartesian base-vectors be ($\mathbf {e}_x,\mathbf {e}_y,\mathbf {e}_z$) in space $\mathbb {R}^3$, where position vector is $\mathbf {r}=x\mathbf {e}_x+y\mathbf {e}_y+z\mathbf {e}_z$ and global coordinates are marked as $r_\text {glob}=(x,y,z)$. Consider a regular, compact surface $\Omega$ in $\mathbb {R}^3$, which has a continuously differentiable parametrization with parameters $p$ and $q$ as

Let us place the origin of the local coordinate system at point $\mathbf {o}(p,q$) on the surface $\Omega$, where two base vectors $\mathbf {f}_1$ and $\mathbf {f}_2$ define tangential plane at the point, and the third base vector $\mathbf {f}_3$ is a local surface normal vector, see Fig. 1. These vectors will be nonzero due to regularity of the surfase. When normalized, these vectors form an orthonormal base ($\mathbf {f}_1,\mathbf {f}_2,\mathbf {f}_3$) and are defined as

The final local base ($\mathbf {e}_1,\mathbf {e}_2,\mathbf {e}_3$) is obtained from the base ($\mathbf {f}_1,\mathbf {f}_2,\mathbf {f}_3$), when the local polarization of the electric field is chosen. Let us set

 

Fig. 1. Vector relations between source point base-vectors on an arbitrary surface and a global origin.

Download Full Size | PDF

Polarization on the surface $\Omega$ can also have externally sourced initial conditions, for example, from a polarization of an incident beam, that creates the electric field distribution on the surface $\Omega$. Let $\mathbf {p}=\mathbf {p}(p,q)$ be a vector that is perpendicular to the plane of local field polarization. Then, local direction $\mathbf {e}_1$ is an intersection of the plane normal to $\mathbf {p}$, and the tangential plane of surface $\Omega$ spanned by ($\mathbf {f}_1,\mathbf {f}_2$). Thus $\mathbf {e}_1$, the base vector along the local electric field, and $\mathbf {e}_2$, the base vector along the local magnetic field, are obtained as

The vectors $(\mathbf {e}_1,\mathbf {e}_2,\mathbf {e}_3$) form a right-handed orthonormal base, which together with the origin $\mathbf {o}$ defines a local coordinate system. Let vector $\mathbf {r}$ be presented in this base as $\mathbf {r}=\bar {x}\mathbf {e}_1+\bar {y}\mathbf {e}_2+\bar {z}\mathbf {e}_3$. Its coordinate triple is marked briefly as $r_{\text {loc}}=(\bar {x},\bar {y},\bar {z})$. When the vector $\mathbf {r}$ is expressed with both $(\mathbf {e}_x,\mathbf {e}_y,\mathbf {e}_z$) and $(\mathbf {e}_1,\mathbf {e}_2,\mathbf {e}_3$) bases, as it is known, the dependence of the coordinates on each other is determined by an orthogonal transformation matrix

Let us denote briefly $\mathbf {\bar {r}}=\mathbf {r}-\mathbf {o}$. Thus, we have presentations of the location $\mathbf {\bar {r}}$ in both local and global coordinate systems. The dependence of the corresponding coordinates can be written using the transformation matrix Eq. (5) as

Let us divide the segment of curved surface $\Omega$ containing the electric field distribution $E_0(p,q)$ into a finite partition of small separate sets $\cup _c\Omega _c=\Omega$. Select the point $\mathbf {o}(p,q$) on differential element $\Omega _c$, the vectors associated with it are $\mathbf {e}_1$ and $\mathbf {e}_2$. Area $\Omega _c$ is projected to the plane spanned by the vectors $\mathbf {e}_1$ and $\mathbf {e}_2$. The resulting area at the plane is marked as $\Omega _t$. Likewise, the function $E_0$ is projected into that plane, limited to set $\Omega _t$ and zero elsewhere. This geometry is presented in Fig. 2. The resulting projection is marked as $E_0^t$ (source point amplitude and phase). In the local coordinate system it holds $\bar {z}=0$ on the piece $\Omega _t$, because this is located on the $\bar {x}\bar {y}-$plane. $E_0^t$ is valid when $\bar {x}\approx 0$ and $\bar {y}\approx 0$.

 

Fig. 2. Surface segment $\Omega _c$ projected into a plane $\Omega _t$.

Download Full Size | PDF

The angular spectrum theory is locally applied on the local coordinate system to function $E_0^t$ valid on piece $\Omega _t$. Let’s observe the field it creates at point $\mathbf {\bar {r}}$, where the coordinates are as in Eq. (6). At a fixed point $\mathbf {r}$ holds [1]

As noted at the beginning of this section, this limit also works physically when considering smooth surfaces. It is even more precise when the wavelength decreases compared to the local radius of the surface. Finally, a further examination is made to obtain a closed form for clause Eq. (10); we present proof in the next section.

Heuristically, when $|\Omega _t|$ shrinks, the function $E_0^t/[E_0^t(\bar {0},\bar {0})|\Omega _t|]$, normed in volume, approaches a Dirac delta $\delta _{(0,0)}$ in the sense of the distribution theory; $E_0^t/[E_0^t(\bar {0},\bar {0})|\Omega _t|]\rightarrow \delta _{(0,0)}$ thus $E_0^t\rightarrow E_0^t(\bar {0},\bar {0})|\Omega _t|\delta _{(0,0)}$. As known, the Fourier transform of Dirac delta is 1; thus, by small $|\Omega _t|$

An electromagnetic beam in the plane wave spectrum representation has polarization vector components along the propagation direction. To account for polarization let’s expand $\mathbf {E}_1(\mathbf {r};\mathbf {o})=E_1(\mathbf {r};\mathbf {o})\mathbf {e}_1$ along plane $(\mathbf {e}_1,\mathbf {e}_3)$ by a component $\mathbf {E}_3(\mathbf {r};\mathbf {o})=E_3(\mathbf {r};\mathbf {o})\mathbf {e}_3$. Let $\mathbf {k}=k_{\bar {x}}\mathbf {e}_1+k_{\bar {y}}\mathbf {e}_e+k_{\bar {z}}\mathbf {e}_3$ be the wave vector, where $k_z=[k^2-k_{\bar {x}}^2-k_{\bar {y}}^2]^{1/2}$. Every single plane wave component in the integral representation of the final incident electric field, like Eq. (14), has to be orthogonal to $\mathbf {k}$. As in [27] we obtain

Hence, a total incident electric field is

Furthermore, the polarization can be defined as

There are two critical observations regarding this derivation. First, $E_1(\mathbf {r};\mathbf {o})\mathbf {e}_1$ is analogous to $-\nabla \times (\mathbf {M}G)$ where the scalar Green’s function ($G$) is replaced by its angular spectrum given by Weyl identity, and where $\mathbf {M}=-2\mathbf {n}\times \mathbf {e}_2$ is a small magnetic dipole [28]. Thus, Eq. (13) is a particular case of the Stratton-Chu equation for the electric field in which only the magnetic sources $\mathbf {M}=-\mathbf {n}\times \mathbf {E}_{surf}$ are presented and the electric ones $\mathbf {J}=\mathbf {n}\times \mathbf {H}_{surf}$ have been removed [29]. Secondly, if the surface $\Omega$ is planar, Eq. (17) and Eq. (19) return at the half-space $z\geq 0$ to the original two-dimensional angular spectrum method.

2.2 Proof of curved boundary integral method

Theorem 2.1 Let the function $E_0$ be continuous in $\Omega$ and let $\mathbf {r}\in \mathbb {R}^3$ be such that $(\mathbf {r}-\mathbf {o})\cdot (\frac {\partial \mathbf {o}}{\partial p}\times \frac {\partial \mathbf {o}}{\partial q})\neq 0$ for all $\mathbf {o}\in \Omega$. The field $\mathbf {E}_1$ is created by Eq. (13) and (14) at point $\mathbf {r}$.

Proof. Define $\epsilon >0$. Let us show first that when the partition $\Omega _t$ is sufficiently small, in other words, the areas $|\Omega _t|$ are small enough regardless of t, by replacing the function $E_0^t$ supported by each projected plane $\Omega _t$ with a constant $E_0^t(\bar {0},\bar {0})=E_0(p,q)$, the error made in Eq. (8) is smaller than $\epsilon /4$. Function $E_0$ is uniformly continuous on the compact segment of surface $\Omega$, and based on this, with sufficiently small $|\Omega _t|$ it holds

Due the monotone of integral, $I(|\bar {z}(p,q)|)\leq I$ on all $t$ regardless of partitions. The upper limit for the error due $E_0(p,q)$ on the Eq. (8) is obtained as

In the above approximation, a constant $E_0(p,q)$ is used on a region of $\Omega$; as a function, it is zero on the plane outside of this region. Let us denote this function simply as $E_0(p,q)$. Next, consider the Fourier transform of the function $E_0(p,q)$. Because function $e^{i|\bar {z}|\sqrt {k^2-k_{\bar {x}}^2-k_{\bar {y}}^2}}$ is integrable, there can be found an origin centered closed disk $\bar {B}$ such that

Let us now consider the error ($err_2$) made in integral Eq. (8) when $\mathcal {F}\left \{E_0(p,q)\right \}$ is replaced by the constant $E_0(p_t,q_t)|\Omega _t|$. Under the estimates Eq. (26) and Eq. (27) it follows

For all $(k_{\bar {x}},k_{\bar {y}})\in \bar {B}$ and ($\bar {x},\bar {y})\in \Omega _t$, we can estimate by continuity of the function

Replacing by $\mathcal {F}\left \{E_0(p,q)\right \}$ with the constant $E_0(p,q)|\Omega _t|$ in Eq. (8) we make an error

Thus, when the partition $\Omega _t$ is dense enough and $\mathbf {r}\in \mathbb {R}^3$ be such that $(\mathbf {r}-\mathbf {o})\cdot (\frac {\partial \mathbf {o}}{\partial p}\times \frac {\partial \mathbf {o}}{\partial q})\neq 0$ for all $\mathbf {o}\in \Omega$, we get by Eq. (25) and Eq. (32)

Finally, by continuity of the integrand, the integral in Eq. (13) exists, and, because of the estimate Eq. (33), for the points $\mathbf {r}\in \mathbb {R}^3$, as in Theorem 2.1, there also exists

The proof above also holds for the component $\mathbf {E}_3(\mathbf {r})=\frac {1}{4\pi ^2}\iint _\Omega \mathbf {e}_3(\mathbf {o})E_0(\mathbf {o}) E_3(\mathbf {r};\mathbf {o})d\Omega$ of the field in Eq. (17); only a few adaptions are needed. Thus, we can also obtain $\mathbf {E}_3$ as a limit of polyhedra approximation and the planar theory.

2.3 Magnetic field representation

Let’s fix $\mathbf {0}\in \Omega$ and consider a point $r$. which has non-zero local coordinate $\bar {z}$. Then we can take a curl of $\mathbf {E}(\mathbf {r};\mathbf {o})=E_0(E_1\mathbf {e}_1+E_3\mathbf {e}_3)$ at $r$. As well known, a magnetic field at $r$ is

Now the real propagation direction of Poynting vector $\mathbf {S}=0.5\text {Re}\left \{\mathbf {E}\times \mathbf {H^*}\right \}$ in half-space $z<0$ can be chosen with the use of sign($\bar {z}$) in the $H_1$ and $H_2$ field components. With the sign($\bar {z}$) in $H_1$ and $H_2$ components, electromagnetic energy propagates away from the surface $\Omega$ on both $\pm \bar {z}$ half-spaces. Without the sign($\bar {z}$) in $H_1$ and $H_2$ fields, the electromagnetic energy propagates through the surface $\Omega$, similarly as it propagates through the focal plane in 2D ASM. If $\bar {z}\neq 0$, we obtain a total magnetic field at $\mathbf {r}$ as

Simliar to the case of the electric field, we can also obtain the magnetic field as a limit of polyhedral approximation of $\Omega$ using the same project tangential plane theory on their faces at point $r$, where $\bar {z}\neq 0$ for all $\mathbf {o}\in \Omega$. The proof of Theorem 2.1 is also valid here; only a few adaptions are needed.

It can be proven that Eq. (17) and (37) for the electric and magnetic fields at all point $\mathbf {r}\in \mathbb {R}^3\setminus \Omega$ satisfies a light condition related to $\Omega$. For example, on sphere $S^2$, all points satisfy that condition. This situation is confirmed in the simulations presented in following sections. The required extended proof is outside the scope of this work and will be presented in a separate article. Its difficulty lies primarily in that integrals do not exist in the ordinary sense when $\bar {z}=0$.

2.4 Repropagation between surfaces

Consider the electromagnetic field $(\mathbf {E},\mathbf {H})$ at some point $\mathbf {g}$ on the repropagation surface $G$. At that point, the direction of the field’s propagation is along Poynting vector $\mathbf {S}(g)=0.5\text {Re}\left \{\mathbf {E(g)}\times \mathbf {H(g)^*}\right \}$. The "differential plane" T, from which the local incident field is thought to be repropagated back, is perpendicular to the direction $\mathbf {S}(\mathbf {g})$, and $\mathbf {g}\in T$. Let’s create a local coordinate system, which origin is fixed on $\mathbf {g}$ and an orthonormal base $(\mathbf {f}_1,\mathbf {f}_2,\mathbf {f}_3)$, where $\mathbf {f}_3=-\mathbf {S}(g) \|\mathbf {S}(g)\|^{-1}$ and $\mathbf {f}_1$ and $\mathbf {f}_2$ are in the differential plane T. Let’s find the direction of $\mathbf {f}_1$ as a real part of projection of $\mathbf {E}$ into the plane T as

Hence, the normalized base vector for the electric field is

The complex vector $\mathbf {E}$ is presented in this base as

Using Eq. (43), we obtain the local coordinates at point $r\in \mathbb {R}^3$ as

Hence, the repropagated electric field is

The repropagated magnetic field can be obtained as a curl of $\mathbf {F}(\mathbf {r})$ as in the subsection 2.3.

3. Results

3.1 Electric and magnetic fields from a source point

The total electric and magnetic fields from surface $\Omega$ are obtained as a superposition of fields from source points. As shown in derivation, when the surface area of the differential source plane approaches zero, the result Eq. (13) is analogous to the angular spectrum presentation of $\mathbf {E}_1(\mathbf {r};\mathbf {o})=-\nabla \times (\mathbf {M}G)$, where $G$ is a scalar Green’s function and magnetic dipole is defined as $\mathbf {M}=-2\mathbf {e}_3\times \mathbf {e}_1$. We expanded the $\mathbf {E}_1(\mathbf {r};\mathbf {o})$ presentation to include the $\mathbf {E}_3(\mathbf {r};\mathbf {o})$ component, which enables propagating polarization state and allows us to compute the magnetic field with via the curl.

Next, we present the radiation patterns $\mathbf {E}(\mathbf {r};\mathbf {o})$ and $\mathbf {H}(\mathbf {r};\mathbf {o})$ from a source point, where the magnetic dipole is positioned to the origin along the y-axis. Fields are integrated over a disk $\bar {B}(\mathbf {0};10k)=\big \{(k_{\bar {x}},k_{\bar {y}}) \ | \ \sqrt {k_{\bar {x}}^2+k_{\bar {y}}^2}\leq \,10 k \big \}$ to obtain propagating and evanescence fields for near-field presentation. Also, the fields are integrated over $\bar {B}(\mathbf {0};k)=\big \{(k_{\bar {x}},k_{\bar {y}}) \ | \ k_{\bar {x}}^2+k_{\bar {y}}^2\leq k^2 \big \}$ to obtain only propagating fields for far-field presentation. See the simulation arrangement in Fig. 3. The electric field polarization is now along the xz-plane, and the field propagates to the positive z-direction. Figures 4–7 illustrate the electric and magnetic field radiation with $\bar {B}=(\mathbf {0},10k)$ and $\bar {B}=(\mathbf {0},k)$ integration limits. As seen from the simulations in Fig. 4–7) the integration limit should be selected according to the application. When the radius of curvature, of the surface on which the distribution is defined, decreases compared to the wavelength or the fields are computed close to the surface distribution, evanescent fields should be accounted for by integrating over a sufficiently wide range in the real $k$-domain.

 

Fig. 3. Source point simulation arrangement, where magnetic dipole $\mathbf {M}$ is positioned along y-axis. The electric field is computed on xy- and xz-planes. The magnetic field is computed in xy- and yz-planes.

Download Full Size | PDF

 

Fig. 4. Electric field $|E_1|$ and $|E_3|$ components, where $|E_2|=0$, and total field $|\mathbf {E}|=|\sqrt {E_1^2+E_3^2}|$ with integration limits $\bar {B}=(\mathbf {0},10k)$ including evanescence fields.

Download Full Size | PDF

 

Fig. 5. Electric field $|E_1|$, $|E_3|$ components, where $|E_2|=0$, and total field $|\mathbf {E}|=|\sqrt {E_1^2+E_3^2}|$ with only propagating fields by integration limits $\bar {B}=(\mathbf {0},k)$.

Download Full Size | PDF

 

Fig. 6. Magnetic field $|H_1|$, $|H_2|$ and $|H_3|$ components at xy-plane, and $|H_2|$ and $|H_3|$ components, with total field $|\mathbf {H}|$ at yz-plane, where $|H_1|=0$. Integration limit $\bar {B}=(\mathbf {0},10k)$ including evanescence waves.

Download Full Size | PDF

 

Fig. 7. Magnetic field $|H_1|$, $|H_2|$ and $|H_3|$ components at xy-plane, and $|H_2|$ and $|H_3|$ components, with total field $|\mathbf {H}|$ at yz-plane, where $|H_1|=0$. Integration limit $\bar {B}=(\mathbf {0},k)$ including only propagating waves.

Download Full Size | PDF

3.2 Field synthesis with an increased radius of curvature

Let’s examine the error created by the curved surface synthesis by comparing the results with physical optics simulations, which utilize electric and magnetic surface currents. These simulations synthesize the field from an electric field distribution on the spherical surface, polarized along the $\mathbf {e}_\theta -$component (standard spherical coordinates). The source distribution $A_0$ is selected to be equal amplitude over $22.5^\circ$ solid angle with sharp edges (tophat) for creating strong interference patterns. The fields are synthesized from a sphere with the fixed radius $\alpha$ and evaluated on a transverse plane fixed with the size of $12.8\times 12.8 \ \alpha$ with $200\times 200$ datapoints. This plane is positioned $3.8 \ \alpha$ distance from the sphere’s origin. Then the ratio between the wavelength and radius of the sphere size parameter $k\alpha =2\pi \alpha /\lambda$ varies from $4\pi$ to $200\pi$. The arrangement and synthesized field with the size parameter $k\alpha =20\pi$ is illustrated in Fig. 8.

 

Fig. 8. The simulation arrangement and global Cartesian electric and magnetic field components of the incident field synthesized from a tophat distribution with $22.5^\circ$ solid angle. In this simulation, fields are computed on a $12.8\times 12.8 \ \alpha$ plane, located $3.8 \ \alpha$ from the origin of the sphere with a size parameter $k\alpha =20\pi$.

Download Full Size | PDF

The difference between the CBIM and physical optics simulations are compared as an average of difference on each datapoint on the evaluation plane as

 

Fig. 9. The amplitude and phase difference between CBIM and physical optics simulations from a curved surface with RoC of $\alpha$. The size parameter $k\alpha$ is decreased from $200\pi$ to $4\pi$.

Download Full Size | PDF

3.3 Repropagation example

To further demonstrate utility, we synthesize electromagnetic fields from a tophat electric field distribution $A_0$, which is polarized along spherical $\mathbf {e}_\theta$ base vector. The $A_0$ is positioned on the $\alpha$ radius sphere’s surface and propagated onto a $10\times 10\alpha$ evaluation plane at the $5\alpha$ distance from the sphere’s origin, see Fig. 10. Then, from these electric and magnetic fields on the evaluation plane, we resynthesize the field back to the spherical surface and compare it to the original field distribution. The repropagation is simulated using $20\pi$ and $200\pi$ size parameters with $\lambda /10$ spacing between the source points. Figures 11–12 illustrate the repropagated field distributions on the spherical surface. The original tophat field distribution $A_0$ has a $22.5^\circ$ solid angle, and the spherical evaluation subsurface has a $30^\circ$ solid angle. In an ideal case, repropagated fields should be distributed only on the original polarization $\mathbf {e}_\theta$ components. Nevertheless, we have a strong tophat distribution on the $\mathbf {e}_\theta$ component with a flat phase and mild fluctuation in the amplitude, which has a sharper tapering with a shorter wavelength. Also, the "leaking" to other spherical components decreases with the shorter wavelength. In many applications, the precise control of the phase is the most crucial attribute.

 

Fig. 10. Electromagnetic beam synthesized from the tophat electric field distribution positioned on the $\alpha$ radius sphere using $40\pi$ size parameter. The synthesized field is evaluated at the $10\times 10\alpha$ transverse plane positioned at $5\alpha$ distance from the sphere’s origin.

Download Full Size | PDF

 

Fig. 11. Amplitudes and phases of the repropagated spherical field components on the spherical surface using $k\alpha =20\pi$ size parameter in the simulations.

Download Full Size | PDF

 

Fig. 12. Amplitudes and phases of the repropagated spherical field components on the spherical surface using $k\alpha =200\pi$ size parameter in the simulations.

Download Full Size | PDF

4. Conclusions

This article introduced the curved boundary integral method (CBIM) for electromagnetic field synthesis based on the angular spectrum method. The CBIM leverages the concept of discretizing an arbitrary surface into infinitesimally small elements and models the radiated field from each of them. Simulations suggest accurate approximation of diffractive fields from surfaces with RoC larger than 2.5 times the wavelength. The CBIM’s advantages include non-singularity at source points, enabling synthesis from any compactly supported and continuous electric field distribution on a regular surface inside the computational domain.

In addition, the four essential points of CBIM are that (1) CBIM is strictly a source method radiating away at both sides from the field distribution, unlike the source-free angular spectrum method. (2) CBIM is analogous to the magnetic dipole presentation of the Stratton-Chu method with polarization perpendicular to the local propagation direction. We added a polarization component to the propagation direction, which enables us to obtain magnetic field presentation via curl. Both polarization components satisfy Maxwell’s equations. (3) CBIM returns to the nominal 2D ASM presentation in the half-space of positive z when the source surface is planar. (4) Each source point radiates a field to other points, modifying the total source distribution. This filters non-physical distribution $E_0$ to a physical one, which always satisfies Maxwell’s equations. On the other hand, the final form of the more complex source distributions should be verified by simulations.

Furthermore, the article explored the repropagation of fields between surfaces, presenting a method using a local coordinate system and transformation matrix. This approach provides a comprehensive framework for simulating field interactions across surfaces.

In conclusion, the development of the curved boundary integral method extends the capabilities of the angular spectrum method by addressing its limitations related to curved surfaces. It offers new possibilities for applications in optics, imaging, and beam-shaping scenarios, where precise control and manipulation of electromagnetic beams are essential. Further research and applications of this method will advance our ability to model and simulate complex wave phenomena.