Surface Flux Transport Modeling Using Physics-informed Neural Networks

SFT modeling is a technique employed to characterize the movement of radial magnetic flux on the solar surface, driven by various physical processes. The temporal evolution of the radial component of the magnetic field on the solar surface, B_r(λ, ϕ, t), where λ, ϕ, and t represent the latitude, longitude, and time, respectively, accounts for various plasma flows on the solar surface as primary drivers of magnetic flux transport. This includes the meridional flow directed poleward, the differential rotation, and the diffusive cancellation arising from flux merging as indicated in Equation (1). Additionally, the emergence of new ARs over time t is incorporated into the model through the inclusion of a source term denoted as S(ϕ, λ, t).

where λ is the latitude (−90° to +90°), ϕ is the longitude (0°–360°), and R_⊙ represents the solar radius. Here, u(λ) represents the advection velocity associated with the meridional circulation and is given by

where the value of u₀ = 12.5 ms⁻¹ and λ₀ = 75° is the latitude above which the meridional flow becomes negligible. The differential rotation Ω(λ) is employed by the latitude-dependent profile (H. B. Snodgrass 1983)

We have chosen diffusion coefficient η = 500 km² s⁻¹ (C. R. DeVore et al. 1984) and the decay constant τ = 5 yr in our calculations (A. Lemerle et al. 2015; K. Petrovay & M. Talafha 2019). The source (S) of the magnetic field differs depending on the situation and will be described in the context of each specific case. The 1D SFT equation obtained by taking the longitude (ϕ) average of Equation (1 on both sides takes the form:

The above equation has been studied for predicting the strength and timing of sunspot cycle (I. Baumann et al. 2006; R. H. Cameron et al. 2012; H. Iijima et al. 2017; L. A. Upton & D. H. Hathaway 2018; P. Bhowmik & D. Nandy 2018; G. Hazra 2021; J. Jiang et al. 2022; E. M. Golubeva et al. 2023).

The SFT equation is a linear differential equation and hence the solutions can be expressed as linear superpositions. This implies that if we have solutions for two or more different initial conditions, the solution for the combined initial conditions can be calculated by superposing the individual solutions. This ensures a basis set of solutions that can be superposed to address any complex scenario. This property can be exploited using Green’s function approach, which is a powerful method for solving linear differential equations. This interesting connection will be discussed in more detail later.

The temporal evolution of the magnetic field using the 1D SFT equation (Equation (4)) is achieved using the Runge–Kutta 2 (RK2) time stepping. The RK2 coefficients K₁ and K₂ for the ith grid point are calculated as

where t_n denotes the time of the nth time step, denotes the magnetic field at the ith grid point for the nth step and Δt being the time step. The temporal update of the magnetic field for the (n + 1)th time step is calculated as

This method ensures the second-order accuracy in time.

The advection part (A_i) is calculated using the upwind scheme with a Van Leer flux limiter (B.V. Leer 1977). The advection update is calculated using the advection flux that depends on the advection velocity (u; see Equation 2) and magnetic field (B) as shown below:

Here, Δλ denotes the latitude grid spacing. The flux values on the left and right grid walls of the ith grid cell are denoted by and , respectively. In order to ensure stability, (Δt) is selected so as to satisfy the Courant–Friedrichs–Lewy (CFL) condition.

The advection flux follows an upwind selection rule:

where . The left and the right values of the magnetic field and are calculated as

Here the harmonic mean slope limiter δ B_i is given by

where .

The diffusion term (R_i) at the ith grid point for the nth time step is constructed as shown below

where η and τ represent the diffusion coefficient and decay rate as described in Section 1.1. To calculate the evolution of B due to the diffusion term (R_i) in the 1D SFT equation, we have used the tri-diagonal solver. The tri-diagonal matrix inversions are performed through the Thomas algorithm (W. H. Press et al. 1992). This method is second-order accurate in space throughout the computational domain.

We have employed the Neumann boundary condition with a zero value at the latitude boundaries, ensuring unimpeded plasma outflow through these boundaries. The selection of the boundary condition is determined by the physical model, which suggests that the magnetic flux will be propelled toward the poles by meridional flows, eventually accumulating and playing a role in polarity inversion during the subsequent cycle. In the case of the 2D SFT equation (Equation (1)), we implemented periodic boundary conditions along the longitudinal axis.

RK-IMEX is a second-order technique (in space and time) that evolves the PDE’s advection term using an upwind Godunov-type explicit approach and solves the diffusion part implicitly (L. Pareschi & G. Russo 2005; S. Kundu et al. 2021). The advection term A_i is calculated using the advection flux as described in the Explicit scheme (see Equations (9)–(13)). The diffusion term R_i is calculated as described in Equation (14). Similar to the Explicit case, Δt is selected using the CFL condition. To evolve by time step Δt to get , we use the following algorithm,

where

The first step is an implicit update of the diffusion term using a tri-diagonal solver. In the second step, the advection flux is explicitly calculated using the previously obtained state , while the diffusion term is treated implicitly. In the final step, the state of the next time step is calculated explicitly using the states already obtained and . The boundary conditions are similar to that mentioned in the Explicit scheme (Section 2.1).

To advance the evolution of the magnetic field in both λ and ϕ directions using the 2D SFT equation (Equation (1)), we employ the dimension-splitting technique (E. F. Toro 1999; D. Gąsiorowski 2012) in both schemes (Explicit and RK-IMEX). This allows the use of the 1D algorithms described above to solve the 2D equation where the fundamental equation can be expressed as a linear sum of operators in each direction independently. The 2D SFT equation (Equation (1)) can be written in the form

where L^λ and L^ϕ represent the operations with λ and ϕ, respectively. The numerical update is performed using dimensional splitting and the Strang operator (G. Strang 1968), as shown below:

In the first step, the magnetic field is updated for L^λ over a duration of to acquire . In the subsequent step, L^ϕ is updated for a full time step of Δt, to obtain . Finally, to complete one full time step update, the magnetic field is once again updated for L^λ over a duration of . The Strang method is employed to ensure a second-order accurate temporal evolution even when extending the individual 1D second-order algorithms to 2D calculations.

Figure 1 illustrates the general architecture of PINNs, which includes the input layer (pink color), hidden layers (red), and output layer (cyan) connected through weights (black solid arrows) and biases (M. Raissi et al. 2019). The input layer comprises the latitude (λ), longitude (ϕ), and time (t), while the output predicts the magnetic field (B) at the specified (λ, ϕ) coordinates at the given time t. Hidden layers consist of fully connected neural networks, forming the overall architecture of the PINN model. Hyperbolic tangent (tanh) is used as the activation function f(z). Here, z represents the weighted sum of values from the preceding neural network layer along with the bias term (I. Goodfellow et al. 2016). The training starts by randomly setting the weights and biases for each connection in the network. For every training point, the network computes the output magnetic field B using PINNs. The loss function is then computed based on satisfying the governing equation (SFT equation), along with initial and boundary conditions. The weights are further updated through backpropagation to minimize the value of the loss function. This iterative process continues until the loss function stabilizes.

Zoom In Zoom Out Reset image size

In PINNs, the computation of derivatives relies on the automatic differentiation method utilizing the network’s structure (weights, biases, and activation function). The outputs are linked to the inputs through activation functions, which have known derivatives. We can determine the derivative of the output (B) with respect to the inputs (t, ϕ, and λ) by applying the chain rule. This makes the PINN method independent of computational grids and allows derivative calculations at any point in the domain. This meshless characteristic of the PINN method is a major advantage, effectively addressing the concerns related to resolution dependence commonly being observed in numerical schemes described in Sections 2.1 and 2.2 above. The absence of any mesh significantly reduces the artificial-intelligence-based algorithm’s vulnerability to numerical diffusion (Q. Hou et al. 2023).

The hyperparameter set (Θ) that includes the network depth and width, number of training points, optimizer, learning rates, and iterations is optimized as described in Section 3.2. In this work, we have employed both the ADAM optimizer (D. P. Kingma & J. Ba 2014) and the L-BFGS optimizer (D. C. Liu & J. Nocedal 1989). ADAM is applied first, followed by L-BFGS. This sequential approach takes advantage of each optimizer's strengths: ADAM provides rapid initial convergence, while L-BFGS ensures precise final optimization. For a given set Θ, the loss function used for the learning process is calculated as:

where f(λ, ϕ, t; Θ) is the residual of the SFT equation given by

Here B_pinn(λ, ϕ, t; Θ) is the predicted value from the PINNs at the point in the domain (λ, ϕ, t). N_ic, N_bc, N_pde, respectively, stand for the number of points selected for the training in the initial, boundary condition, and domain. B_init(λ, ϕ) denotes the initial condition used in the simulation, whereas represents boundary conditions at boundary points λ_bc and ϕ_bc. Since the problem has two spatial dimensions, the boundary conditions for both λ and ϕ have been incorporated separately and subsequently added to get the total loss function for the boundary condition. w_ic, w_bc, w_pde, respectively, are the weights given to each term in the loss function during the training process. Each term in the loss can have a different order of magnitude, which results in unequal training of individual terms. The inclusion of the weight in the loss function forces the network to converge to a unique solution following the given constraints. The mean square error (MSE) has been used as the loss function. Since the 1D SFT equations (Equation (4)) do not include terms in ϕ, the corresponding input and the loss term () are omitted while simulating the 1D case.

Hyperparameter optimization plays a crucial role in assessing the SFT equation using PINNs. The overall optimization of the code for both 1D and 2D has been performed using Gaussian process (GP)-based Bayesian optimization. In this method, an initial set of hyperparameters is assumed to begin the optimization process. The network is optimized to calculate the test and train loss using these values. Once the losses are calculated, a GP-based regressor is employed to predict the next set of hyperparameters using the current state. Hence there is a two-level optimization to obtain the best set of hyperparameters (P. Escapil-Inchauspé & G. A. Ruz 2023). A brief description of the optimization technique used is given in the Appendix. We optimized all the hyperparameters including the number of network layers, the number of nodes in a layer, the number of points used for training (initial, domain, boundaries), the iterations required in each optimizer, the learning rate, the activation function, and the weights for the loss function (domain, boundaries, and initial). A hyperparameter list for both 1D and 2D simulations obtained after the optimization is given in Table 1. It is worth noting that in the case of 1D simulation, the SFT parameters described in the introduction part are used along with the source described in the 1D section later. For the 2D case, we take a BMR as described in case 2D. We tested the code for different initial conditions and found that the hyperparameter set is consistent. The same has been verified using multiple simulations with different initial positions, orientations, and polarity reversed.

Table 1. Hyperparameters Used in 1D and 2D Simulations

Parameters		1D		2D
Hidden layers	[nodes] × layers	[41] × 10		[32] × 8
Number of training points	Initial (Nic)	2787		28838
	Boundary (Nbc)	2356		20022
	Domain/PDE (Npde)	87460		173567
Optimizer		ADAM	L-BFGS	ADAM	L-BFGS
Weights on the loss	Initial (wic)	0.60	0.32	0.17	0.25
	Boundary λ ()	0.01	0.04	0.16	0.25
	Boundary ϕ ()	-	-	0.45	0.25
	Domain/ PDE (wpde)	0.39	0.64	0.22	0.25
Learning rate		0.0022	-	0.00052	-
Iterations		35016	14166	199495	7904

Note. This table describes the hyperparameter set including the hidden layer configuration, number of training points for different phases and optimizer choices. The values are obtained using the GP-based Bayesian optimization. In both cases, training is conducted in two steps: first with ADAM, and then with L-BFGS. Weights assigned to loss components, learning rates, and iterations for both optimization steps are also mentioned.

Download table as:  ASCIITypeset image

To comprehend the relationship between the size of the training data set and the model, several training sessions, each with a varying number of training points, have been carried out. The L1 error for both the test and train sets, corresponding to various sizes of the training data set for the 2D case, is illustrated in Figure. 2. Despite initial fluctuations, the errors stabilize after ∼2 × 10⁵ training points for both test and training data sets. This observation serves as a confirmation that the model optimization is successful and that the provided architecture consistently yields minimal errors.

Zoom In Zoom Out Reset image size

To validate the numerical (RK-IMEX and Explicit) and the PINN methods, we use the analytical solution prescribed by C. R. DeVore et al. (1984). The analytical method takes the meridional flow and diffusion terms in the 1D SFT equation into account (first and second terms on the right-hand side of Equation (4)). The relative effect of meridional flow and diffusion is represented by . We use β = 10 for our validation following the range of values adopted by Leighton (R. B. Leighton 1964; C. R. DeVore et al. 1984).

The model assumes an initial state B(μ, t = 0) as shown below

where and θ = 90°–λ. A special advection profile u(θ) is considered for this specific case given by

Here u₀ = 12.5 m s⁻¹ denotes the amplitude of the velocity profile. The diffusion coefficient η = 500 km² s⁻¹ is used as described in Section 1.1. This gives an exponentially decaying analytical solution as given below

with decay time constant .

In Figure 3, we compare the results obtained in the PINNs and RK-IMEX with the analytical solution. The absolute error between the PINNs and analytical solution is in the range (0.2%–1.3%). As expected for any numerical scheme, the RK-IMEX method shows higher accuracy in the analytical solution on increasing the number of grid points (Ngrid; see Figure 3(b)). This is further supported by the error plot, where the L1 norm demonstrates a decreasing trend while increasing Ngrid from 256 to 2048 bins. PINNs exhibit the lowest error even compared to the highest Ngrid simulation using the RK-IMEX scheme (Figure 3(a)). Because of its mesh independence, the PINN solution exhibits fewer fluctuations in the error.

Zoom In Zoom Out Reset image size

The total flux is calculated as

The comparison of flux calculated for the southern hemisphere is plotted in Figure 3(c). It is observed that the flux loss decreases on increasing Ngrid. The error in flux from the PINN solution is ∼1.07%, which is lower than that of the RK-IMEX method with the highest Ngrid (see Figure 3(d)). The PINN solution produces the best match for the flux compared to the numerical method. This shows that the flux loss in the numerical method is higher than in the PINN method. The flux loss could critically affect the solutions produced using SFT simulations, especially when the aim is to calculate the polar flux at the end of a solar cycle.

To verify the stability of the chosen hyperparameters, we performed a sensitivity test (ST). This involved randomly altering some of the hyperparameters and retraining the network. Table 2 illustrates the perturbed hyperparameters and their resulting magnetic and flux errors for various sets. The change in hyperparameters is expressed in percentages, where an increase (decrease) is indicated by ↑ (↓). The magnetic field error was found to be high in all the sets in comparison to that of the optimized hyperparameter set (error = 0.61%) (Figure 4(a)). We found that the flux errors range from 1% to 14% as depicted in Figure 4(b). The flux errors are lower in some of the sets (ST1,4,6,8 and 10), although corresponding errors in the magnetic field are higher (as compared to the error for the optimized set of hyperparameters). This exercise indicates the consistency and stability of the selected hyperparameter set.

Zoom In Zoom Out Reset image size

To validate the numerical (RK-IMEX and Explicit) and PINN methods, we evolve the BMRs after initializing at various latitudes (0° and 12°) as shown in Figures 5(a) and (d). We utilized a resolution of Δλ = Δϕ = 1° for our simulations. Once the BMR is initialized, the average magnetic field is calculated along the longitude to make the 1D magnetic profile (Figures 5 (b) and (e)). The initial magnetic field distribution in latitude (λ) is fitted using the function given below

Zoom In Zoom Out Reset image size

The fitted parameters are A = 14.88, B = 0.14, C = −3.65, D = −7.07, λ_c = −2.33, and σ_c = 3.51 for the first case (0°), and A = –14.88, B = 0.14, C = −2.24, D = –7.07, λ_c =−0.67, and σ_c = 3.51 for the second case (12°). This distribution is then evolved using the 1D SFT equation without the source term. The BMR and parameters are chosen from A. R. Yeates (2020). The evolved magnetic field, along with the respective initial states and BMRs are plotted in Figures 5(c) and (f).

When the BMR is located at the equator (0°), each magnetic polarity is carried away in opposite directions, resulting in distinct polarities at each pole. The flux cancellation is very low and the flux loss is mainly due to the diffusive and decay losses. When the BMR is initiated in the northern hemisphere (12°), the magnetic field mostly gets pushed toward the corresponding pole. The magnetic field intensity gets reduced due to the flux cancellations between different polarities advecting together. There is a very small amount of field that has crossed the equatorial region because of the initial diffusion. The solutions from the RK-IMEX and PINN schemes obtained in the present study are compared with the standard numerical scheme from A. R. Yeates et al. (2023); hereafter referred to as Yeates scheme (see Figure 6 in their paper). We observe that the features obtained in the evolution of BMR <B_r (t) > are similar to columns 3 and 4 of Figure 5 of this work. We have further estimated the L1 norm obtained between the solution from both the schemes (RK-IMEX and PINN) and that of the Yeates scheme. For cases 1 and 2, the difference is 9% and 14%, respectively, between the RK-IMEX and Yeates schemes. Such a deviation is expected as the Yeates scheme is first-order accurate, whereas the RK-IMEX scheme has second-order accuracy. The difference between the PINN method and the Yeates scheme is 11% and 15% for cases 1 and 2, respectively. Additionally, the difference between the RK-IMEX and PINNs method is found to be ∼6% in both cases 1 and 2. This is expected because of the lesser numerical diffusion in the PINN method compared to the other two numerical methods.

Zoom In Zoom Out Reset image size

We evolved the magnetic field using the 1D SFT equation (Equation (4)) by including the meridional flow, diffusion, decay, and source terms using both numerical schemes (Explicit and RK-IMEX) and PINNs. The parameters used for advection, diffusion, and decay are described in Section 1.1. The source term S(λ, t) represents a smooth distribution that characterizes the probability of preceding and following polarities emerging on the solar surface. This representation is akin to the approach employed by M. Dikpati et al. (2006) and consists of a pair of rings featuring opposite magnetic polarities. The path of the new ARs followed by a complete cycle (Λ₀) is determined by the quadratic approximation established based on extensive observations of solar cycles (J. Jiang et al. 2011).

where P = 11 yr is the time period of a typical solar cycle. The source term is given by the equation

where k = ± 1 is used alternating between even and odd cycles and A_m = 0.003 is the amplitude of the source term (K. Petrovay & M. Talafha 2019). The time profile S(t) of solar activity in a typical cycle was determined from the average of many cycles as (D. H. Hathaway et al. 1994)

with a = 0.00185, b = 48.7, c = 0.71, and t_c representing the time since the last cycle minima. The latitudinal profile S₀(λ; Λ₀, δ λ) is a Gaussian given by

where the full width at half maximum δ λ = 6°. As the solar cycle progresses, sunspots tend to appear at higher latitudes during the cycle’s peak and migrate toward lower latitudes as the cycle declines. This phenomenon contributes to the latitudinal separation (Δλ) observed in the rings and is a result of Joy’s law as shown below.

We defined our simulation box to span from −70° to +70°, excluding the polar regions, to avoid the singularity associated with magnetic flux and term in the SFT equation. The simulation was carried out for 24 yr, with the initial condition set as . We have employed boundary conditions as mentioned in the Explicit scheme (Section 2.1).

The comparison of the solution of the PINNs and both numerical methods (Explicit and RK-IMEX) with the same initial condition and empirical equations are given in Figure 6. The solution obtained using the RK-IMEX (Ngrid=2048) and PINN schemes are plotted in Figures 6(a) and (b), respectively. We have also obtained solutions for Ngrid values 128 and 512 using RK-IMEX. Figure 6(c) illustrates the outcomes achieved by PINNs and RK-IMEX for different values of Ngrid at a particular time instance (t = 6.384 yr). This temporal point is indicated by black dashed lines in Figures 6(a) and (b). For the sake of brevity, we have not plotted the solution obtained from the Explicit scheme as it is similar to the RK-IMEX solution. The results demonstrate a noticeable alignment between the outcomes obtained through the RK-IMEX and PINN method. However, it is important to recognize that differences exist in the solutions, particularly near the boundary as evident from the inset plot in Figure 6(c). There is a distinct shift in the time of polarity inversion near the poles (see black rectangular boxes in Figures 6 (a) and (b)). The L1 norm of both the numerical schemes (RK-IMEX and Explicit) with respect to the PINN solution is shown in Figure 6(d). By increasing the number of grid points, the numerical solutions for both the RK-IMEX and Explicit schemes gradually converge toward the PINN solution. Both the RK-IMEX and Explicit schemes exhibit a similar order of accuracy, whereas RK-IMEX has a lower error value.

The loss function for training the PINNs for 1D is shown in Figure 7. The weighted loss plots the loss value for which the optimization was performed that includes the weights in each term (see Equation (22) and Table 1). The unweighted loss is calculated after dividing out the weights from the weighted loss. The ADAM optimizer effectively reduces the loss initially by considering the initial magnetic field and the source terms that are given higher weights. The L-BFGS optimizer further minimized the loss of each term to achieve improved accuracy. The loss function tends to stabilize and saturate as the training progresses toward the later iterations, indicating that the optimization process has converged to a satisfactory solution. This convergence assures that the results obtained are reliable and reflect a good level of accuracy in our model.

Zoom In Zoom Out Reset image size

We use the numerical (Explicit and RK-IMEX) and PINN methods to evolve a BMR planted at a point on the θ–ϕ plane using the 2D SFT equation (Equation (1)). The model parameters including the advection profile (u), diffusion coefficient (η), and decay constant (τ) are chosen as described in Section 1.1. The initial condition is taken as a BMR with positive and negative circular regions for the respective polarity AR emergence. We follow the prescription by A. A. van Ballegooijen et al. (1998) to artificially introduce this BMR as described in Equations (36)–(39).

Here, β_± are the heliocentric angles between each point (θ = 90°–λ, ϕ) and the center of each polarity region. ρ₀ represents the angular separation between the center of the positive and negative polarities (θ_±, ϕ_±).

An angular width of σ = m ρ₀ is assumed around the center of each polarity. For small β_±, the BMR is similar to a Gaussian distribution with this angular width to account for the contribution of the respective polarities. A linear sum of the contribution from each polarity is used as the initial magnetic field distribution. To illustrate, we present results using m = 0.4 and (θ₀, ϕ₀) = (75°, 270°) as parameters for generating the artificial AR. The evolution of this initial state is performed using Equation (1) without including the source term for 110 days. We implemented periodic boundary conditions along the longitudinal axis and outflow conditions along the latitudinal boundaries as described in Section 2.1.

The comparison between the solution obtained from PINNs and the RK-IMEX (Δλ = Δϕ = 05) scheme is depicted in Figure 8. The BMR generated using the above parameters is plotted in Figure 8(a). The evolved states for the numerical and PINN methods are shown in Figures 8(b) and (d), respectively. The PINN result exhibits a high degree of similarity to the RK-IMEX solution, indicating strong agreement between the two (Figure 8(c)). The mean relative error between the solutions is ∼1% (with mean absolute error < 0.0012 G) for the resolution Δλ = Δϕ = 05. The error calculated with respect to the PINN solution is observed to fall on increasing the resolution, similar to the 1D case. Hence, it can be inferred that the numerical solution on increasing the number of grid points approaches the PINN solution.

Zoom In Zoom Out Reset image size

To assess the additivity of the magnetic field with respect to the SFT equation, we introduced two BMRs at different locations on the θ–ϕ plane. This approach is possible because the SFT equation is a linear differential equation (see Section 1.1). The initial seed BMRs (BMR1 and BMR2) are centered at (θ₁, ϕ₁) = (75°, 270°) and (θ₂, ϕ₂) = (90°, 234°), respectively (Figure 9(a)). Initially, these BMRs are evolved independently using the PINN method (Figures 9(c) and (e)). Figure 9(d) demonstrates the sum of the individually evolved states. The PINN method is also used to solve the SFT equation for both the BMRs initiated together. The solution obtained by evolving the magnetic field simultaneously is shown in Figure 9(b). These two solutions are then compared to understand the additivity of magnetic field for multiple BMRs. The difference between the solutions is ∼1% (Figure 9(f)). This demonstrates that the evolution of multiple BMRs can be effectively accomplished using the PINNs approach by individually evolving each BMR and then summing their contributions. The PINN method is also tested for initial conditions with reversed polarity and has been demonstrated to work effectively.

Zoom In Zoom Out Reset image size

The applicability of the models to observed magnetogram data is assessed by evaluating them using data from SOHO-MDI. SOHO is a joint project between the European Space Agency (ESA) and NASA. We have used the “fd_M_96m_01d” data series that consists of the full disk magnetogram of (1024 × 1024) grids with (4″) resolution in a spectral range of 6767.8 Å ± 190 mÅ and a cadence of 96 minutes (P. H. Scherrer et al. 1995). AR 7978 is used as an illustration of the application and can be replicated in any other AR. This AR was observed in the southern hemisphere during the solar minima between cycles 22 and 23 (1996). It persisted for multiple CRs (1911–1916), allowing for the examination of its development for at least two rotations after the peak of flux emergence (A. Ortiz et al. 2004; L. van Driel-Gesztelyi & L. M. Green 2015). Moreover, it has distinct north and south polarities, making it an ideal candidate for studying the applicability of the PINN model. To generate the magnetic field of the initial condition, we applied the BMR prescription as described in Equations (36)–(39). The PINN and RK-IMEX methods are used to evolve the BMR over time (see Figure 10). To visually analyze the alignment of PINN’s BMR solution and the actual data, the BMR region is marked as contours in the data. The contour traces the magnetic field flow over the data, verifying the evolution. The total unsigned flux is calculated for both the PINN and RK-IMEX solutions with respect to the data. The flux error for PINNs ranges from ∼1% to 10%, whereas it is ∼3%–11% for the RK-IMEX method. The flux error is lower for the PINN method as we obtained in the 2D case. The exact values of the flux error may vary based on the parameters of the SFT modeling and the algorithm used for generating the initial BMR.

Zoom In Zoom Out Reset image size