• Home
• Sitemap
Article View

## Technical Paper

Applied Science and Convergence Technology 2019; 28(3): 66-78

Published online May 31, 2019

https://doi.org/10.5757/ASCT.2019.28.3.66

## Visualizing Electromagnetic Vector Fields in Matter using MATHEMATICA

Yong-Dae Choi, and Hee-Joong Yun*

Department of Microbial & Nanomaterials, Mokwon University, Daejeon 35349, Republic of Korea College of Techno-Sciences, Mokwon University, Daejeon 35349, Republic of Korea

Correspondence to:*E-mail: heejy@mokwon.ac.kr

Received: December 28, 2018; Accepted: February 27, 2019

This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-CommercialLicense (http://creativecommons.org/licenses/by-nc/3.0) which permits unrestricted non-commercial use, distribution,and reproduction in any medium, provided the original work is properly cited.

The quantities of dynamic physics including gravitational field and electromagnetic field are accurately described using vector differential operators. Electromagnetism is notably more conceptual in matter compared to vacuum because abstract, three-dimensional fields that are sometimes difficult to visualize are involved in most analysis. A visual representation of abstract vector fields in matter is invaluable to students or researchers working in this field and may also assist teachers in teaching electromagnetism to physics or engineering students. We successfully visualized the most fundamental concepts of the electromagnetic vector calculated and estimated using vector differential operators in Mathematica. This visualization based on vector calculations can be used as a starting platform for the further exploration of electromagnetic vector fields.

Keywords: Electromagnetic vector field, Mathematica, Vector differential operators, Mathematica simulation, Vector field visualization

Electromagnetic vector fields are composed of elaborate mathematical structures and abstract concepts described as vector differential operators [15]. It is often conceptually challenging for physicists to intuitively understand the associated physics phenomena and distinguish between them. To learn about the changes in the electric displacement $D⃗$ in a medium, including polarization $P⃗$ and electric quadrupole Qα, β, or magnetic intensity $H⃗$ in the medium including the magnetization vector $M⃗$ in the field $B⃗$0 [5,6], the vector fields must be calculated and visualized to understand their characteristic features. According to our educational experience, the physicists must introduce a visualizing portfolio in curriculum to examine and test these concepts of physics. Several tools are available for the visualization of electromagnetic field vector such as: MATLAB’s interactive tools [7,8], Maple’s Classroom Tips [9], and Mathematica demonstrations [1014]; that are effectively used to plot vector fields [79,12]. However, they are limited to certain fields in vacuum [9,13] or special materials [1518], and do not specifying the constitutive relations in matter. The electric and the magnetic fields $E⃗$ and $B⃗$ respectively, are determined every where in space if, all sources, charge ρ and current density $J⃗$, are specified. For a small number of definite sources, the determination of the fields is a tractable problem. However, for most cases of macroscopic aggregates of matter, the $D⃗$ and $H⃗$ are determined by the macroscopic sources $P⃗$, $H⃗$, Qα, β and similar moment densities of the material medium in the presence of an applied field. The derived fields $D⃗$ and $H⃗$ are expressed in terms of $E⃗$ and $B⃗$ using the constitutive relations in the material medium [5,18]. Maxwell equations do not preclude the possibility that one or both of quantities ɛ, μ, be negative. When μ > 0, Faraday’s law in vector differential form, $k⃗$ × $E⃗$ = μω$H⃗$ denote a right-handed orthogonal triad whereas μ < 0, {$k⃗$, $E⃗$, $H⃗$} constitute a left-handed triad, which is explained by the vector differential mode exactly. We presented the capability of Mathematica for vector analysis and visualizing vector fields in previous works [19,20]. In this work, we have presented a platform visualizing vector fields added by the field densities using Mathematica. The platform will manipulate visualizing of five vector fields: $E⃗$ = −$∇⃗$ϕ, $D⃗$ = ɛ$E⃗$, $B⃗$ = $∇⃗$×$A⃗$, $H⃗$ = −$∇⃗$ϕ*, and $B⃗$ = μ$H⃗$ in matter. It is composite of several modules to visualize the electromagnetic vector fields and specify the electromagnetic characteristics of the comparable modes of the material. Therefore, it is potentially useful as a pedagogical tool for the electromagnetism in the higher class of undergraduate. It is assumed that the reader has some familiarity with vector differential operators in electrostatics. The basic background of electrostatics (EM) is reviewed briefly in Section II. In Section III, the plots of vector fields are produced in Mathematica and demonstrated using graphics. In Section IV, a summary is provided based on the platform of vector fields. The reader who is conversant with vector differential mode in EM may skip Section II. In Appendix, few segments of Mathematica 11.0 code are presented to aid the reader’s comprehension of vector fields plotting.

### 1. Electric vector field from scalar potential

(1) Electric field in vacuum

Electrostatics is the study of the time-independent distributions of electric charges and electromagnetic fields in electrodynamics [16]. If an arbitrary charge q is placed in an electric field $E⃗$, it will experience a force $F⃗$ given by $F⃗$ = q$E⃗$. The electric field exerted to test charge q can also be expressed as an integral over the charge density instead of the total charge: Q = ∫ρ($r′→$)dv′

$F→(r→)=14πɛ0qQ|r→-r′→|2e^ r=q14πɛ0∫V(r→-r→′)V|r→-r→′|3ρ(r→′)dv′︸E→(r→)=-q∇→(14πɛ0∫V1|r→-r→′|ρ(r→′)dv′)︸ϕ(r→)▶ E→=-∇→ϕ(r→), with ϕ(r→)=14πɛ0∫Vρ(r′→)V|r→-r′→|dv′,$

where dv′ is used to designate an element of volume in the charge distribution and V denote the entire volume occupied by the charge distribution. We get Gauss’s law in vector differential operator from Eq. (1) by the divergence theorem;

$∇→·E→=ρɛ0 or ∇→ϕ·n^=-σɛ0$

As since $E⃗$ = −$∇⃗$ϕ, Gauss’s law guide to Poisson equation $∇2ϕ=-ρɛ0$. In regions of no charge, Poissons’s equation reduce to Laplace equation: ∇2ϕ = 0. The general solution of Laplace equation in spherical coordinates for a problem with azimuthal symmetry is as follows [5,6]:

$▶ ϕ(r,θ)=∑l=0∞[Alrl+Blr-(l+1)]Pl(cosθ),$

where Pl(cosθ) represents the Legendre polynomials. Coefficients Al and Bl can be determined from the boundary conditions.

(2) Constitutive relation in dielectrics: $D⃗$ = ɛ0$E⃗$+ $P⃗$

For the potential of an arbitrary localized charge distribution ($r′→$) in power of $1τ$, we can write potential at point $r→$13,5] as :

$ϕ(r→)=14πɛ0∫(1τ)ρdv′=14πɛ0∫[1r∫ρdv′+1r2∫r′ cosθρdv′+1r3∫(r′)2(23cos2θ-12)ρdv′+…]$

Here we use the binomial expansion of Legendre polynomials: $1τ=1r∑n=0∞(r′r)nPn(cosθ)$. The first term is called the monopole term, 2nd the dipole term, and 3rd quadrupole term. In this case we take the dipole term as a continuous charge distribution. The potential due to the dipole moment can be written succinctly as:

$ϕp(r→)=14πɛ01r2r^.∫r→′ρdv′=14πɛ0p→·r→r3,$

where $P⃗$ = $r′→$ρdv′ is called the dipole moment of the distribution. A convenient measure of this effect is the polarization $P→-limΔv→01Δv∑ip→i$, then $P⃗$ = P⃗dv [1,2]. And apply the vector identity $∇→·P→(1r)=1r∇→·p→+p→·∇→(1r)$. Then the potential in the dielectrics with $P⃗$ instead of $P⃗$ in Eq. (5):

$ϕ(r→)=14πɛ0∫vP→·r→r3dv=14πɛ0∫vP→·∇→(1r)dv=14πɛ0[1r∫SP→·n^da-∫v1r(∇→·P→)dv]=14πɛ0∫S1rσbda+14πɛ0∫V1rρbdv,$

wheres σb$P⃗$ · is the bounded surface charge and ρb = −$∇⃗$·$P⃗$ is the bounded volume charge due to the polarization of the dielectric media. When applying Gauss’s law to a region containing the charges embedded in a dielectric, the polarization charges (ρb) as well as the charge (ρ) must be included. Gauss’s law in dielectric medium is

$∇→·E→=1ɛ0(ρ+ρb)=1ɛ0(ρ-∇→·P→)⇒ ∇→·(ɛ0E→+P→)=∇→·D→=ρ,$

where $D⃗$ = ɛ0$E⃗$ + $P⃗$ is a new macroscopic field vector, the electric displacements $D⃗$, that act as the constitutive relation in dielectrics. $D⃗$ having the same unit as $P⃗$ charge per unit area (C/m2), for $∮D⃗$·n̂da = Q.

### 2. Magnetic vector field from the potentials

(1) Magnetic vector field: $B⃗$ = $∇⃗$ × $A⃗$

Using the Biot-Savart law, the force $F⃗$2, exerted on circuit 2 of current I2 owing to the influence of circuit 1 in current I1 is as follows:

$F2→=μ04πI1I2∮1∮2dl→2×[dl→1×(r2→-r1→]|r2→-r1→|3=I2∮2dl→2×μ04πI1∮1dl→1×(r→2-r→1)|r→2-r→1|3=I2∮2dl→2×B→(r2)$

Next, using the vector identity $∇⃗$ × (&phi$F⃗$) = $∇⃗$ϕ×$F⃗$ + ϕ$∇⃗$ × $F⃗$, we can write as follows:

$B→(r2)=μ04π∫v1J→(r1→)×(r→2-r→1)|r→2-r→1|3dv1=-μ04π∫v1J→(r→1)×∇2(1|r→2-r→1|)dv1=∇→2×μ04π∫v1J→(r→1)|r→2-r→1|dv1=∇→2×A→(r2→)$

where the calculation $∇⃗$2J(r1) = 0. It confirms that $B⃗$(r2) = $∇⃗$2 × $A⃗$(r2) with $A⃗$(r2)

$▶ A→(r→2)=μ04π∫v1J→(r→1)|r→2-r→1|dv1.$

The vector potential $A⃗$($r→2$) may be applied to current circuits by making the substitution $J⃗$dv = $I⃗$dr and use an expand in powers of r1/r2 to get $|r2-r1|-1=1r2(1+r1·r2r22+…)$ and vector identity $dr→1(r→1·r→2)=12(r→1×dr→1)×r→2$ [1]. Then Eq. (10) is

$A→(r→2)=μ0I4π∮dr→1|r→2-r→1|=μ0I4π1r2∮dr1+1r23∮dr→1(r→1·r2→)+…≃μ04π[12∮(r→1×dr1→)×r2→r23]=μ04πm→×r→2r23,$

where $m=12∮(r→1×dr→1)$ is the magnetic dipole moment of the circuit.

(2) Constitutive relation in magnetic media: $B⃗$ = μ0 ($H⃗$ + $M⃗$)

Now we define a macroscopic vector point function magnetization: $M→=limΔv→01Δv∑imi$. Each volume element ΔV of magnetized matter is characterized by a magnetic dipole moment: $M⃗$ = M⃗dv′ to explore the influence of matter on the magnetic field. We can write $B⃗$ with $ϕM*(r→)$ and $M⃗$ from Eq. [11]:

$B→=∇→×A→=μ04π∫v0∇→×(M→(r→′)×(r→-r→′)|r→-r→′|3)dv′=μ04π∫v0[M→(r→′)∇→·(r→-r→′)|r→-r→′|3-(M→(r→′)·∇→)(r→-r→′)|r→-r→′|3]dv′=μ0M→(r→)-μ0∇→(14π∫v0M→(r→′)·(r→-r→′)|r→-r→′|3dv′)︸ϕM*(r)B→(r→)=μ0M→(r→)-μ0∇→ϕM*(r→)=μ0(M→(r→)+H→(r→)),$

where magnetic scalar potential $ϕM(r)*$ is

$ϕM(r)*=14∫v0M→(r→′)·(r→-r→′)|r→-r→′|3dv′)$

Here an auxiliary magnetic vector field, $H⃗$, is introduced as $H→=1μ0B→-M→$, i.e., $B⃗$ = μ0($H⃗$ + $M⃗$), this is the constitutive equation in magnetic medium calculated from the magnetic scalar potentials : $H→(r→)=-∇→ϕM*(r→)$. Using the vector identity $∇→′·(M|r→-r→′|)=∇→′(1|r→-r→′|)·M→+1|r→-r→′|∇→′·M→$, we can present $ϕM*(r→)$ with the pole densities σM and ρM [13,5].

$ϕM*(r→) =14π∫V0M→(r→′)·∇→′(1r→-r′→)dv′ =14π∫V0∇→′·(M→(r→′)|r→-r→′|)-1|r→-r→′|∇→′·M→(r→′)dv′ϕM*(r→) =14π∫S0σM(r→′)|r→-r→′|da′+14π∫V0ρM(r→′)|r→-r→′|dv′,$

where σM($r′→$) = $M⃗$($r′→$) is the surface density of magnetic pole strength and ρM ($r′→$) = −$∇⃗$′·$M⃗$($r′→$) is the magnetic pole density obtained from an integral over the entire volume of material V0 in Eq. (14), which contrast with the potential from the dipole moment ϕ($r′→$) in Eq. (6).

(3) Maxwell equation in matter

The Maxwell’s equations describing the vector fields in macroscopic media can be written in the vector differential form at the MKSA unit system [13,5] as :

$∇→·D→=ρ, ∇→×H→=J→+∂D→∂t$$∇→·B→=0, ∇→×E→+∂B→∂t=0$

where $D⃗$ = ɛ0$E⃗$ + $P⃗$, $H→=1μ0B→-M→$. In dielectric media, polarization $P⃗$ = ɛ0χ$E⃗$ is defined and that $D⃗$ = ɛ0(1+ χ) $E⃗$ = ɛ$E⃗$, where χ is the electric susceptibility of a material. The dielectric constant K is written by $K=ɛɛ0$ and $D⃗$ = ɛ$E⃗$ = 0$E⃗$. In magnetic material, magnetization is described by $M⃗$ = χm$H⃗$, then $B⃗$ = μ0(1+ χm) = $H⃗$ = μ0Km$H⃗$ = μ$H⃗$ where $Km=1+χm=μμ0$. Km is the relative permeability. As $∇⃗$·$B⃗$ = 0, and $∇2ϕM*=0$, if μ is at least piecewise constant in each region, the magnetic scalar potential $ϕM*$ also satisfies Laplace’s equation [5,6]. To visualize vector fields, we must determine the potentials from the Eqs. (1), (6), (10), and (14) with zonal harmonics. Eq. (3) is applicable in spherical boundaries and azimuthal symmetry applications. Thus we can calculate and plot the vector fields appropriately on the medium environments.

### 1. Plotting E→ = −$∇⃗$ϕ

For the visualization of vector fields, Mathematica supports many procedures including VectorFieldPlot, GradientFieldPlot, and ContourPlot. VectorFieldPlot [{fx, fy}] plots a vector function {fx, fy}, and GradientFieldPlot[(f)]plots the gradient of scalar function f. VectorFieldPlot[{x,y}] and GradientFieldPlot[(x2 + y2)] plots are shown in Fig. 1(a) column. They are simpler but instructive examples that show the same vector fields results from PlotVectorField[{x,y}] and Plot GradientField[(x2 + y2)]. $E⃗$ = −$∇⃗$ϕ is shown in bottom row of Fig. 1. We have used the userdefined procedure heejyPlot (in Mathematica code 3) to plot and manipulate the electric and magnetic fields with the prescribed scalar potential (In[451–458] in Mathematica code 2) appropriate to the material medium more effectively. The procedure includes the modules GradientFieldPlot and ContourPlot, and material parameters available to the prescribed potentials are assigned. When we assigned a potential the procedure produce a vector field platform as in Fig. 2. The plots of −$∇⃗$ϕ are produced with the GradientFieldPlot function in Mathematica. Complete modules are available in Ref. [24].

### 2. Plotting $B⃗$ = $∇⃗$ × $A⃗$

It is helpful to calculate vectors in both Cartesian and cylindrical coordinates in Mathematica as in Mathematica code 1. If circuit current is confined in wire the vector potential $A⃗$2($r→$) is determined by Eq. (10), then we can calculate magnetic induction $B⃗$ with curl of vector potential $B⃗$($r→$) = $∇⃗$ × $A⃗$z($r→$) and plot vector fields in both 2D and 3D graphics with VectorFieldPlot3D in Mathematica as shown Fig. 3. The calculated $B⃗$ field is $BX∞={-2I0y(x2+y2),2I0x(x2+y2),0}$, in Cartesian coordinates while $BC∞=[0,-2I0r,0]$

in cylindrical coordinates in Mathematica code 2. This indicates $B⃗$ = $∇⃗$ × $A⃗$ fields have x and y component in Cartesian while only θ component in cylindrical coordinates when $A⃗$z($r→$) have a z component, which is the characteristics of $B⃗$ = $∇⃗$ × $A⃗$. If, instead of being confined to wire, the current exists in a medium, then $A⃗$z($r→2$) is represented with the magnetic dipole moment $m→=12∮(r→1×dr→2)$ as in Eq. (11).

### 3. Prescribed potentials for vector fields

(1) Polarized sphere with polarization $P⃗$ in vacuum

The electric fields generated by a uniformly polarized sphere of radius a with polarization $P→=limΔv→01Δv∑ip→i$ or $p⃗$ = $P⃗$dv. The potential due to $P⃗$ is given by Eq. (6),

$ϕ(r→,θ)=14πɛ0∫s1rσbda=14πɛ01r∫sPcosθda,$

where σb = $P⃗$ · = Pcosθ, and the volume bounded charge density ρb = −$∇⃗$ · $P⃗$ = 0 because $P⃗$ is the constant polarization vector. Laplace’s equation in spherical coordinates is solved by zonal harmonics with spherical boundaries and azimuthal symmetry as follows [1,2];

$ϕ1(r,θ)=∑l=0∞[Blr-(l+1)]Pl(cosθ),r≥a$$ϕ2(r,θ)=∑l=0∞[AlrlPl(cosθ)], r≤a$

At r = a, the boundary conditions: ϕ1(r, θ) = ϕ2(r, θ) give the relation $A1a1=B11a1+1$, and $∇→ϕ1r-∇→ϕ2r=-σ(θ)ɛ0$ by Eq.(2). From these we can get $A1=P3ɛ0$ by the “eyeball” method; σ(θ)b = $P⃗$ · = PCos [θ]. Thus, the potentials as follows:

$▶ ϕp1(r,θ)=P3ɛ0a3r2Cos [θ], r≥a$$▶ ϕp2(r,θ)=P3ɛ0rCos [θ], r≤a$

Next, we obtain the electric fields inside and outside of the uniformly polarized sphere; E1 = −$∇⃗$ϕ1(r, θ) E2 = −$∇⃗$ϕ2(r, θ), and plot as Fig. 5(a).

(2) Grounded spherical conductor in electric field $E⃗$0

As a canonical problem in physics, we considered a conducting sphere radius a held at constant (let zero) potential and placed in a uniform field $E⃗$0 [1,5,6,10]. The potential can be written as an expansion with Legendre polynomials by Eq. (3). The coefficients {Al, Bl} were obtained from the boundary conditions. The electric field should be $E0→$ as r –> ∞. So that the potential must be $ϕ1=-∫E0→·d→l=-E0r Cos [θ]$. Then $ϕ(r,θ)=∑l=0∞AlrlPl(cosθ)=-E0r Cos [θ]$, therefor all Al = 0 except A1 = −E0. At r = a ϕ = 0, Alal + Bla−(l+1) = 0⇔B1 = −A1a3 = E0a3. Thus we determined the potential:

$▶ ϕc1(r,θ)=-E0(rCos [θ]-a3r2Cos [θ])$

Next, we get $E⃗$1 = −$∇⃗$ϕ(r, θ) and plot as Fig. 5(b).

The field intensity is determind using huefunc $[Norm [epvx [z,y]][epvxmax]]$ with ColorFunctions. The field intensity platform is manipulated with the materials (K values) and medium environment (initial fields) and calculations of the vector field intensity at the vector points (r = 0.75a,1.15a) are shown in Fig. 4.

(3) Dielectric sphere in an initially uniform electric field $E⃗$0

The potential of the dielectric sphere with a radius of a in initially uniform electric field $E⃗$0 may be expressed as a sum of zonal harmonics in Eq. (18) [1,2]. As r –> ∞, $E⃗$ –> $E⃗$0 and the potential on surface should be continuous. The potential is presented by two harmonics of the lowest-order for the outside and inside of dielectric sphere, respectively as follows:

$ϕ1(r,θ)=A1rCos [θ]+C1r-2Cos [θ]], r≥a$$ϕ2(r,θ)=A2rCos [θ]+C2r-2Cos [θ]], r≥a$

The constants A1, A2, C1, and C2 determined from the boundary conditions as A1 = −E0, C2 = 0, $A2=-3E0K+2$, and $C1=K-1K+2E0a3$. Then we determine the electric scalar potentials:

$▶ ϕe1(r,θ)=E0(K-1)(2+K)a3r2Cos [θ]-E0rCos [θ], r≥a$$▶ ϕe2(r,θ)=-3E02+KrCos [θ], r≥a$

For $D⃗$ = ɛ$E⃗$ = 0$E⃗$, we obtained the dielectric scalar potentials ϕd1 and ϕd2 as ϕd1(r, θ) = ɛ0ϕe1, ϕd2(r, θ) = 0ϕe2. We get $E⃗$1 and $E⃗$2 from the ϕe1 and ϕe2 in an initially uniform field of $E⃗$0 as follows:

$E→1=-∇→ϕe1={(E0+K-1K+2a3r3Cos [θ],(-E0+K-1K+2a3r3)Sin[θ],0}, r≥a$$E2→=-∇→ϕe2=3E0K+2{Cos [θ],-Sin [θ],0}, r≤a$

We can write $E→2=3K+2E0e^z$ in Cartesian coordinate, where êz = Cos[θ]erSin[θ]êθ. This can be rewriten as $E→2=(1-K-1K-2)E0e^z$, which represent the depolarization (Clausius-Mossotti factor $K-1K+2$) [1,2]. The molecular polarizability is $α=3ɛ0NK-1K+2$, where N is number of atoms or molecules per volume. Thus, characteristics of the $E⃗$2 inside the dielectric give us important quantitative information about the structure of the atoms [21]. The fields $E⃗$ and $D⃗$ inside and outside of the dielectric sphere on the K values are shown in Fig. 5.

(4) Uniformly magnetized sphere in a nonpermeable medium

If the current density vanishes in some finite region of space, $∇⃗$ × $H⃗$ = 0 by Ampere’s law (Eq. (16)). This implies that we can introduce a magnetic scalar potential ϕM such that ∇2ϕM = 0. As a common practical situation, we consider the hard ferromagnets with magnetized $M⃗$. For such material, fixed magnetization $M⃗$($x⃗$). Then, $∇⃗$ × B = μ0$∇⃗$ × ($H⃗$ + $M⃗$) = 0 becomes a magnetostatic Poisson’s equation ∇2ϕM = −4πρM with a magnetic charge density ρM if there is no boundary surface. The magnetic scalar potential ϕM due to the magnetization $M⃗$ is given by [1,5]:

$ϕM(r,θ)=13Ma2rCos [θ]$

Thus, the magnetic potential is

$▶ ϕM1(r,θ)=13Ma3r2Cos [θ], r≥a$$▶ ϕM2(r,θ)=13Mr Cos [θ], r≥a$

$M⃗$ = z is parallel to the z-axis. For the magnetic induction in the magnetized sphere, we may estimate a modified the magnetized scalar potential $ϕM2*$ in the sphere satisfying $B→2*=-∇→ϕM2*=μ0(H→+M→)=μ023M→$, thus, we can determine $ϕM2*$ as follows:

$▶ ϕM2*(r,θ)=-μ023Mr Cos [θ], r≤a$

Then $H→2=-∇→ϕM2=-13M→$ and $B→2*=-∇→ϕM2*=μ023M→$. Here, we note that $B→2*$ is parallel to $M⃗$ while $H2→$ is antiparallel to $M⃗$ as shown in Fig. 6(a).

(5) Uniformly magnetized sphere in uniform magnetic field $B⃗$0

Consider a sphere of linear magnetic material with a radius a and a permeability μ in space containing an initially uniform magnetic field $B⃗$0. The boundary conditions are applied to obtain the parameters of the Eq. (23): (i) as r –> ∞, $B⃗$1 = $B⃗$0 from −μ0$∇⃗$ϕ2 = $B⃗$0, we get $A1=-B0μ0$. (ii) as r –> 0, C2 –> 0 and (iii) $B⃗$1r = $B⃗$2r at r = a, from those ${C1→a3B0(-1+Km)(2+km)μ0,A2→-3B0(2+Km)μ0}$ [1,5]. Thus, we get the magnetic potentials:

$▶ ϕm1(r,θ)=B0μ0[(-1+Km)2+Kma3r2-r]Cos [θ], r≥a$$▶ ϕm2(r,θ)=-B0μ03(2+Km)r Cos [θ], r≤a$

Then, we get $B⃗$1 = μ0$H⃗$1 = μ0(−$∇⃗$ϕm1) and B2 = μH2 = μ(−$∇⃗$ϕm2) in the spherical system as follow:

$B→1=B0{(-1+Km)2a3Cos [θ](2+Km)r3+Cos [θ],(-1+Km) a3Sin [θ](2+Km)r3-Sin [θ],0}, r≤a$$B→2=Km{3B0Cos [θ]2+Km,-3B0Sin [θ]2+Km,0}, r≥a$

where $B→2=Km3B02+Kme^z$ parallel to the pole. Vector fields of $B⃗$1 and $B⃗$2 are plotted on the Km values in the Cartesian coordinate in Fig. 6. Thus, we can calculate $B⃗$1 and $B⃗$2 with magnetic scalar potentials conveniently rather than with vector potentials.

(6) Prescribed magnetic scalar potentials

We can visualize the electromagnetic vector fields with heejyPlot user defined procedure composed of GradientFieldPlot and ContourPlot using the prescribed scalar potentials (▶ marked). We then summarize the electromagnetic scalar potentials and vector field functions in Mathematica code 1 that are utilized in the heejyPlot function. We can use the above prescribed eight potentials (In[451–458]) of Cartesian coordinates to plot the vector fields in the heejyPlot procedure. By changing the parameters (K, Km, M, P) of the medium, plot the vector fields can be plotted repeatedly using the prescribed potentials shown as in Figs. 5 and 6.

### 4. Vector field plotting with scalar potentials

(1) Electric field plotting in the matter

Visualization of the electric fields in matter is presented two material medium: (i) the electric fields of the near the material in vacuum and (ii) in a uniform electric field $E⃗$0. In vacuum, consider of the electric field produced by a uniformly polarized sphere of radius a and polarization $P⃗$. To plot the electric fields ($E⃗$1, $E⃗$2) inside and outside the uniformly polarized sphere in vacuum where no other electric fields are present, we employed the prescribed potential potp12[y_, z_]. For plotting the electric fields ($E⃗$1, $E⃗$2) and displacements ($D⃗$1, $D⃗$2) fields inside and outside the material sphere we have assigned the prescribed potentials: potc1 [y_, z_], pote12 [y_, z_], potd56 [y_, z_] to the heejyPlot procedure. The prescribed potentials are described in the Mathematica code 1. According to the selection of potential and parameter procedure plot the vector fields array as in Fig. 5. The adjustment of table index steps and plotting of additional vector field frames are possible by the index of table as in Mathematica code 2. It is evident that $E→2=-P3ɛ0e^z$ opposite to $P⃗$ and so $E⃗$2 called the depolarization field [22,23]. Given that $D⃗$2 = ɛ0$E⃗$2 + $P⃗$ while $E⃗$2 = −$∇⃗$ϕe2 and $P⃗$ = ɛ0χ$E2⃗$, we can examine $D⃗$ and $E⃗$ outside and insid sphere for the K values shown as Fig. 5. Visualization of the electric and the magnetic vector fields is presented using both vector field plotting and vector field density plotting with vector point field intensity values (see Field intensity at (II) in Figs. 5 and 6). The vector field densities are plotted using $huefunc=Norm [vfunc]vfmax$ with ColorFunction as in Mathematica code 2. These plot effectively reveal the features of the vector fields in matter. We can confirm the features of the vector fields by calculating the vector points inside and outside the sphere in the same platform. The intensity {E1 > E2, E1 > E2, E1 > E2, D1 < D2} is shown in Fig. 5 and the array {H1 > H2, B1 <B1, H2> E2, B1 < B2} is shown in Fig. 6.

(2) Magnetic field plotting in the matter

For the plots of magnetic intensity ($H⃗$1, $H⃗$2) and magnetic induction ($B⃗$1, $B⃗$2) fields in vacuum, we employed the prescribed potentials: poth42 [y_, z_], potb43 [y_, z_] described in Mathematica code 1. To plot the magnetic intensity ($H⃗$1, $H⃗$2) and magnetic induction ($B⃗$1, $B⃗$2) in the unform magnetic field $B⃗$0, the prescribed potentials: potm77 [y_, z_], potm78 [y_, z_] are utilized heejyPlot. The plot frames can be expanded by controlling the table index of the platform similar to the plotting of electric fields. In Fig. 6(a) we see that $H→2=-M3e^z$ is opposite to magnetization vector $M⃗$ = z, that is the demagnetizing field, a result that in accord with $E→2=-P3ɛ0e^Z$ in Fig. 5(a). All the visualization of the vector fields was done by plotting the vector fields instead of the illustrations described in most of texts [15]. It confirms that the constitutive relations: $B⃗$ = μ0($H⃗$ + $M⃗$) and $D⃗$ = ɛ0$E⃗$ + $P⃗$ are graphically represented in Figs. 5 and 6.

(3) Constitutive vector fields platform

Plotting the vector fields in different medium environments, we determined the eight functions for the scalar field potentials as presented in Mathematica code 1. Array of the vector fields is produced by selecting the prescribed potential and proper parameters, K, Km, P, M, E0, B0 etc.. Next, an array of vector fields is automatically produced based on the suggested indexing. In this procedure eight arrays of vector fields are possible. Therefore, it is possible to test and examine the constitutive relations in both vector calculation and graphics mode simultaneously in the same platform. The difference between them could nor shown until vector field plots are produced in different medium environments. We have summarized the vector fields by visual representation of six fundamental vector fields as shown in Fig. 7. The platform presents the plotting of five kinds of vector fields and vector field density: $E⃗$ = −$∇⃗$ϕ, B = $∇⃗$ × $A⃗$, $H⃗$ = −$∇⃗$ϕ*, D = ɛ0$E⃗$+ $P⃗$, and $B⃗$ = μ0($H⃗$ + $M⃗$), those are most fundamental vector fields in electromagnetics. Visualizing vector fields in Mathematica demonstrate the enhanced effectiveness of vector plots (eps, pdf) are scaled up in the graphics. Recently, artificial media i.e., metamaterials exhibit qualitatively new electromagnetic response functions which can not be found in the nature. Metamaterials are characterized by simultaneously negative permittivity and permeability, it referred to left-handed materials (LHMs). Their physical characteristics are also described with Maxwell equation Eq. (16) in vector differential mode [17,18]. Some aspects of Mathematica codes plotting vector fields in figures are presented in Appendix. Complete Mathematica codes plotemvfs.nb plotting the plots in the paper allows users to manipulate vector plots that are available from the site in Ref. [24]. Even if Mathematica is not installed, the user may use a plotemvfs.cdf file instead in CDF Player (free download) [25].

We have visually represented the most fundamental electromagnetic vector fields with the scalar and vector potentials using plotemvfs in Mathemtica. In this work, eight user-defined functions for the prescribed potentials, were assigned, that were directly applied to the procedure heejyPlot to simultaneously plot the vector fields. This enables the examining and testing of the vector fields in different media environment and confirm the constitutive relations: $D⃗$ = ɛ0$E⃗$ + $P⃗$ and $B⃗$ = μ0($H⃗$ + $M⃗$) graphically In this procedure, we manipulated only the vector field plotting of $E⃗$ and $B⃗$ separately unless they interacted with each other. Thus, we only considered the time independent terms in the Maxwell’s equations Eq. (15). The visualization of electromagnetic vector fields was facilitated using the Poynting or Jones vectors in other studies [26,27]. Herein, we focused on visualizing the fundamental electric and magnetic vector fields separately and manipulate the plot of the vector fields in cross-checking with vector differential operators. According to our results, the plotemvfs platform could be a useful foundation for physics students and researchers to conceptualize, visualize, and quantize the law of physics knowledge in the knowledge and information converging process.

The authors would like to thank Professor Taewan Kim for helpful discussion on the visualization of vector field and pedagogical strategies in physics class

The Mathematica codes used for the plots shown in Figs. 46 to interpret the producing the vector fields. However, with these codes the vector fields could not be produced. The vector fields are produced by the complete plotemvfs procedure (Ref. [24]).

### Mathematica code 1

In[431]: = eq4[y_,z_] = phic1, ϕc1 (Eq. (22))In[432]: = eq31[y_,z_] = phip1, ϕp1 (Eq. (20))In[433]: = eq32[y_,z_] = phip2, ϕp2 (Eq. (21))In[434]: = eq1[y_,z_] = phie1, ϕe1 (Eq. (25))In[435]: = eq2[y_,z_] = phie2, ϕe2 (Eq. (26))In[436]: = eq41[y_,z_] = phiM1, ϕM1 (Eq. (30))In[437]: = eq42[y_,z_] = phiM2, ϕM2 (Eq. (31))In[438]: = eq43[y_,z_] = phiM2s, ϕM2* (Eq. (32))In[439]: = eq5[y_,z_] = phid1, ϕd1 (Eq. (25))In[440]: = eq6[y_,z_] = phid2, ϕd2 (Eq. (26))In[441]: = eq70[y_,z_] = phim1, ϕm1 (Eq. (33))In[442]: = eq80[y_,z_] = phim2, ϕm2 (Eq. (34))In[443]: = eq7[y_,z_] = μ0*phim1, ϕm1 (Eq. (33))In[444]: = eq8[y_,z_] = μ*phim2, ϕm2 (Eq. (34))In[451]: = potc1[y_,z_]: = If[y2+z2>a,eq4[y,z],0];In[452]: = potp12[y_,z_]: = If[y2+z2>a,eq32[y,z],eq31[y,z]];In[453]: = pote12[y_,z_]: = If[y2+z2>a,eq2[y,z],eq1[y,z]];In[454]: = poth42[y_,z_]: = If[y2+z2>a,eq42[y,z],eq41[y,z]];In[455]: = potb43[y_,z_]: = If[y2+z2>a,eq43[y,z],eq41[y,z]];In[456]: = potd56[y_,z_]: = If[y2+z2>a,eq6[y,z],eq5[y,z]];In[457]: = potm77[y_,z_]: = If[y2+z2>a,eq80[y,z],eq70[y,z]];In[458]: = potm78[y_,z_]: = If[y2+z2>a,eq8[y,z],eq7[y,z]];

### Mathematica code 2

In[501]: = heejyPlot[potential,xlim,ylim,opts_]  : = Module[{plot1,plot2}, plot1 = GradientFieldPlot[potential  ,xlim,ylim,PlotPoints-> 18  ,DisplayFunction->Identity]; plot2 = ContourPlot[potential  ,xlim, ylim,ContourShading ->False  ,PlotPoints->250]] Show[plot1,plot2,opts]; heejyPlot[potential ,opts_]: = heejyPlot[potential,  ,{x,1.1,−1.1},{y,1.1,−1.1} opts];In[502]: = Table[  kd = {1.0, 2.3,3.5,4.3,5.0,6.1,28.4,80.1 };  kemat = {"K = 1.00 Air","K = 2.3 Polyethylene",  "K = 3.5Nylon","K = 4.3 Quarz","K = 5.0 Glass",  "K = 6.1NaCl","K = 28.4Alcohol","K = 80.1Water"};In[503]: = Ec1 = Show[heejyPlot[potc1[y,z]  ,{z,−2.2,2.2},{y,−2.2,2.2},text1,text62  ,DisplayFunction ->Identity],circle1  ,DisplayFunction ->\$DisplayFunction  , ImageSize -> 161]; Ep2 = Show[heejyPlot[potp12[y,z],ibid.]] Ein = Show[heejyPlot[pote12[y,z],ibid.]] Din = Show[heejyPlot[potd56[y,z],ibid.]]In[510]: = Row[{Ep2,Ec1,Ein,Din}  ,{i,2,5,5}]// FramedIn[520]: = Row[{Ep2,Ec1,Ein,Din}  ,{i,1,8,6}]// FramedIn[523]: = epv1 = -Grad[phip1{r,θ,ϕ},"Spherical"]  /.{r->1.15,P->5.5 10−11} //EvaluateIn[524]: = epv2 = -Grad[phip2{r,θ,ϕ},"Spherical"]  /.{r->0.75,P-> 5.5 10−11}//EvaluateOut[510] = {-Graphics-} (Fig. 4)Out[520] = {-Graphics-} (Fig. 5) Out[523] = {2.2102,0,0}(epv1) Out[524] = {−2.07063,0,0}(epv2)
Fig. 1. (Color online) Plots of the vector fields (top row) and gradient of scalar fields (bottom row). Plot (a) Vector{x,y} is generated using the function VectorFieldPlot[{x,y}] and (x2 + y2) (bottom row) is generated using GradientFieldPlot[(x2 + y2)] in Mathematica. The other plots are gernerated using the same method used in (a).
Fig. 2. (Color online) Plots of E = −φ of the two monopole potentials generated by function GradientFieldPlot in Mathematica.
Fig. 3. (Color online) Plots of vector potential Az(r) and vector field B = × Az(r) in 2D and 3D. These vector fields were visualized using the function VectorFieldPlot3D in Mathematica.
Fig. 4. (Color online) Vector field intensity platform for epvx[z,y] using huefunc[Norm[epvx[z,y]]epvxmax] with ColorFunctions shows the results as epvx[1.15a, 0.75a] = {2.12012, −2.07063} and ecvx[1.15a, 0.75a] = {4.3516, 0}, respectively.
Fig. 5. (Color online) Plots of vector fields (I) and the vector field densities (II) of E and D for a dielectric sphere of radius a. (I) The vector field E[z, y] of a polarized sphere with P = 5.5 × 10−11 C/m2 and E[z, y] and D[z, y] for three dielectric constants, K = {∞, 2.3, 2.3}, in an initially uniform electric field E0 = 2.15 V/m. (II) The vector field densities of E[z, y] and D[z, y] under the same conditions as (I). The vector field densities are estimated from the maximum fields of {Epmax = 1.0 V/m, Ecmax = 3.15 V/m, Edmax = 3.15 V/m, and Ddmax = 3.15/Cm2}.
Fig. 6. (Color online) Plots of the vector fields (I) and the vector field densities (II) of H and B for a magnetized sphere of radius a. (I) The vector field H[z, y] for a magnetized sphere of M0 = 0.97 T, the vector field B[z, y] for a magnetized sphere of M0 = 77,190 T in vacuum, and the vector fields, H[z, y] and B[z, y] for a magnetic material with Km = 2.2 in an initially uniform vector field B0. (II) The vector field densities of H[z, y] and B[z, y] under the same conditions as in (I). The vector field densities are estimated from the maximum fields Hmax = 4π ×105 A/m and Bmax = π T, respectively, to compare the fields in a platform. M0 = μ0M, where M is the saturation magnetization of the ferromagnet [1].
Fig. 7. (Color online) Vector field platform visualization of five vector fields from the vector potentials and scalar potentials: (a) E = −ϕ, (b)B = × A, (c) D2 = ɛ0E2 + P in an initially uniform electric field E0, (d) H = −ϕ* and B2 = μ0(H + M) side by B2*, and norms (d-norms) corresponding to (d) fields. The sphere of M in B2 is embedded in a nonpermeable medium while the sphere of B2* in an initially uniform magnetic field B0. All the plots are selected from the Figs. 1–7. While one scale up the graphics on the platform examines the vector graphics more detail. The number to the left side is the time (seconds) to plot the 15 fields in Mathematica (CPU (dual-core):i3-6100@3.7GHz).
1. JR. Reitz, FJ. Milford, and RW. Christy, Foundation of Electromagnetic Theory, (Addison-wisley, Reading, 1992), pp. 9, 63 97-126, pp.218-282.
2. DJ. Griffiths. Introduction to electrodynamics, (Prentice Hall, Eaglewood Cliffs, 1989), p.118, p.166, p.176, Chaps.2 and 4.
3. DK. Cheng. Fundamentals of Engineering Electromagnetics (Addison-Wesley, San Francisco, 1993).
4. MNO. Sadiku. Elements of Electromagnetics (Oxford University Press, New York, 2001) p. 173.
5. JD. Jackson. Classical Electrodynamics, (JohnWiley & Soms, 1975), p.90, p.196.
6. GB. Arfken, and HJ. Weber, Mathematical Methods for Physics 4th ed., (Academic Press, San Diego, 1995), p. 510 and p. 748
7. O. Ozgun, and L. Sevgi, IEEE Antennas and Propagation Magazine. 57, 113 (2015).
8. HG. Espiosa, and DV. Thiel, IEEE Antennas and Propagation Magazine. 17, 140 (2017).
9. RL. Zimmerman, and FI. Olness. Mathematica for Physics, 2nd ed. (Addison-Wesley, San Francisco, 2002), Chap 8.
10. PT. Tamm. A Physicist’s Guide to Mathematica (Academic Press, San Diego, 1997), Chaps. 4 and 5.
11. Wolfram Research. The Mathematica 11.0 Documenter Center (Wolfram Media, 2017).
12. Wolfram Research. The Mathematica 11.0 Standard Addon Packages (Wolfram Media, 2017).
13. Y. Shibuya, http://demonstrations.wolfram.com , Dielectric Sphere In A Uniform Electric Field (2018).
14. YJ. Kim, Appl Sci Converg Technol. 25, 77 (2016).
15. DH. Park, and NM. Kim, Appl Sci Converg Technol. 25, 154 (2016).
16. V. Veselago, L. Braginsky, V. Shklover, and C. Hafner, J Comput Theo Nanosci. 3, 1 (2006).
17. G. Kim, and B. Lee, J Electromag Engin Sci. 16, 67 (2016).
18. HJ. Yun, New Phys: Sae Mulli. 50, 134 (2005).
19. YD. Choi, and HJ. Yun, J Korean Vacuum Soc. 21, 69 (2012).
20. B. Nathalie, H. Tap, and B. Daniel, Molecular Phys. 102, 783 (2004).
21. NA. Pertseb, and H. Kohlstedt, con-mat. 0603762 (2006).
22. Q. Zhang, and I. Ponomareva, New Phys: Sae Mulli. 15, 043022 (2013).
23. plotemvfsprogram http://www.heejy.site/helicity.html (2019).
24. YD. Choi, and HJ. Yun, J Korean Phys Soc. 67, 792 (2015).
25. TW. Kim, and HJ. Yun, J Korean Phys Soc. 69, 697 (2016).