Optimal control for a coupled spin-polarized current and magnetization system

This paper is devoted to an optimal control problem of a coupled spin drift-diffusion Landau–Lifshitz–Gilbert system describing the interplay of magnetization and spin accumulation in magnetic-nonmagnetic multilayer structures, where the control is given by the electric current density. A variational approach is used to prove the existence of an optimal control. The first-order necessary optimality system for the optimal solution is derived in one space-dimension via Lagrange multiplier method. Numerical examples are reported to validate the theoretical findings.


Introduction
The classical theory of micromagnetism describes the dynamics of ferromagnetic materials occupying a domain in the absence of electric currents and for constant temperature (below the Curie temperature). The process of magnetization is described by so called Landau-Lifshitz-Gilbert (LLG) equation, see [13,17]. A more general approach is to consider the interactions between spin accumulation and the magnetization on magnetic and nonmagnetic multilayer structures [12,24,25,29,30,34]-which has wide applications in various magnetic devices, e.g., in magnetic tunnel junctions and magnetic domain walls [26,31]. Moreover, a number of technological applications of these phenomena have been seen, e.g., in racetrack memories, magnetic vortex oscillators [21,23]. In this paper, we study the optimal control problem for a coupled spin drift-diffusion Landau-Lifshitz-Gilbert system on a magnetic multilayer. A formal description of our problem is as follows.
Let D,D ⊂ R d (d = 1, 2, 3) be two bounded Lipschitz domains such that D ⊂D, and let ∂D and ∂D be their boundaries, respectively. For any T > 0, we denote D T  which is deduced from the Landau-Lifshitz energy E(m). For simplicity, in this work E(m) is the exchange energy • γ 0 and α are the gyromagnetic ratio and the nondimensional empiric Gilbert damping parameter, respectively. • D 0 is the diffusion coefficient. The parameters λ 1 and λ 2 are the characteristic length of the spin-flip relaxation and the mean free path of an electron, respectively. • The parameters β and β are the nondimensional spin polarization parameters of the magnetic layers.
In case of j = 0, and s = 0, (1.1)-(1.3) reduces to the standard LLG equation, the well-posedness of which is studied, e.g., in [6,8,16,20,33] and references therein. In [12], the authors have employed a Galerkin approximation to prove the existence of global weak solution of (1.1)-(1.3) in three space-dimensions. By using energy methods, the authors in [14] established the existence of a global smooth solution of spin-polarized transport equation in two space-dimensions; see also [25] for the existence of a smooth solution in one dimensional case. Recently a decoupled timemarching scheme was analyzed and its unconditional convergence of the integrator towards a weak solution of the underlying problem was established in [3].
An optimal control problem subject to LLG equation has been studied by Dunst et al. in [9]. The authors in [9] have shown the existence of an optimal solution and derived its necessary first-order optimality system in the one space-dimension. Moreover, they have shown convergence of the time semi-discrete optimality system towards the optimality system of the original problem. We also refer to see [4,5] for optimal control type problems subject to LLG equation.
Our goal is to study an associated optimal control problem to (1.1)-(1.3). Let m : [0, T ] × D → S 2 be a given function and let (1.5) for π = (m, s, j) with some nonnegative constants κ, δ 1 , δ 2 and δ 3 . The term (m(T )) represents the terminal payoff. In Definition 2.2 it is given as a Lipschitz continuous function defined on L 2 . In Section 5 for numerical simulations we choose (u) = u −m(T ) 2 L 2 . We find an admissible weak solution π * = (m * , s * , j * ) which minimizes the cost functional, i.e., F (π * ) = min π F (π) subject to (1.1)-(1.3). (1.6) A minimizer π * of (1.6) is constructed via the variational method. For a minimizing sequence of admissible weak solutions π n = (m n , s n , j n ), we employ uniform bounds along with Fatou's lemma to prove existence of a minimizer of (1.6). Once a minimizer π * of (1.6) is obtained, we ask for its necessary optimality systemwhich may be used for a gradient descent method to numerically approximate π * . To deduce necessary first-order optimality system of π * , we need higher regularity of the solutions-which is why, we consider more regular control j in (1.6) and restrict our study to one space-dimension, see, e.g., proof of Lemma 4.2. Due to the presence of nonlinear non-Lipschitz drift in (1.1), the classical Pontryagin's maximum principle is not immediate. We use the Lagrange multiplier method to deduce the necessary first-order optimality system corresponding to minimizer π * in one space-dimension. Due to the coupled system of state equations, it is more challenging to show that π * is a regular point of some appropriate function (cf. Lemma 4.4) -a crucial step for applying the Lagrange multiplier theorem.
The remaining part of this paper is organized as follows. We detail the technical framework, and state the main result in Section 2. Section 3 is devoted to show the existence of an optimal solution of (1.6). In Section 4, we first deduce improved stability properties of the weak solution in one space dimension and then use these estimates to obtain a necessary first-order optimality system for the optimal solution π * via the Lagrange multiplier method. Moreover, improved regularity properties of the adjoint variables and the optimal control are obtained. Computational studies for the switching dynamics are reported in the final section.

Technical framework and statement of the main results
Throughout this paper, we use the letter C to denote a generic positive constant, which may take different values at different occurrences. In the sequel, we denote by iv) m(0, ·) = m 0 and s(0, ·) = s 0 in the trace sense.
Note that it follows from the formulation of (1.1) that |m| is constant. (This is a well-known property of magnetization which states that below the Curie temperature, the magnitude of the magnetization is constant.) Hence we assume |m| = |m 0 | = 1. The existence of global weak solutions for (1.1)-(1.3) is detailed in [3,12]. We assume in this paper that the given data are smooth enough, so that uniqueness holds; see [7,Theorem 1.3]. The following theorem holds: Theorem 2.1 Let s 0 ∈ H 1 (D) and m 0 ∈ H 1 with |m 0 | = 1 a.e. in D, and D 0 :D → R + be a measurable function such that for some positive constants D * and D * . Then for any j ∈ L 2 (H 1 (D)), there exists a weak solution (m, s) of (1.1)- (1.3) in the sense of Definition 2.1. Moreover, the following estimates hold: Thanks to the above theorem, the set U ad of all triplets π = (m, s, j), where j belongs to L 2 (H 2 (D)) and (m, s) is the weak solution to the corresponding problem (1.1)-(1.3), is non-empty. The reason to require more smoothness in j is to obtain more smoothness in the solution (m, s), as will be explained later. With this at hand, we rewrite the minimization problem (1.6) as follows.

Definition 2.2
Let the assumptions of Theorem 2.1 hold true. Letm : [0, T ] × D → S 2 be a given function, and be a given Lipschitz continuous function on L 2 . A tuple π * = (m * , s * , j * ) ∈ U ad is said to be a weak optimal solution of (1.6) if We finish this section by stating the main results of this article, the proofs of which will be presented in Sections 3 and 4.3. The first theorem states the existence of a weak optimal solution π * . Theorem 2.2 There exists at least one weak optimal solution π * ∈ U ad of (1.6) in the sense of Definition 2.2.
In case of one spatial dimension, optimal solution π * may satisfy first-order optimality conditions. For d = 1, we consider the following equation for spin accumulation s :D T → R 3 , see [25], where the spin current J is a vector defined by with the given electric current density j :D T → R. Here, ·, · denotes the Euclidean inner product in R 3 . The next theorem derives the first-order optimality system to be solved to simulate the optimal solution (m * , s * , j * ) numerically. The second-order optimality system is not easy to derive due to lack of regularity of the solution of the control problem.
Thus, thanks to (3.1) and (3.2), there exist subsequences of {s n }, {m n } and {j n } (not relabeled) such that ) and a.e.D T We can achieve the unit length condition for m * in the following way. Let φ 0 ∈ C ∞ 0 (D). Take the test function φ = m * φ 0 in iii) of Definition 2.1. Then one has and therefore |m * (t, x)| = |m 0 (x)| = 1 a.e. D T . Consequently, π * = (m * , s * , j * ) ∈ U ad . It remains to show that π * is a minimum. Observe that F is measurable, nonnegative and lower semi-continuous convex function. Thus, by using Fatou's lemma, we have i.e., F (π * ) = inf π ∈U ad F (π). This completes the proof.
Remark 3.1 One may formulate the control problem for H 1 -control and show the existence of an optimal solution and hence the optimal control. In Theorem 2.2, we have taken H 2 -valued control, which play a crucial role to deduce optimality conditions for optimal solution π * .

Optimality system for d = 1
Theorem 2.2 ensures the existence of an optimal solution (m * , s * , j * ) of (1.6) in any spatial dimension d = 1, 2, 3. In order to deduce the necessary optimality system associated with the tuple (m * , s * , j * ), one needs higher regularity for the solution, and therefore needs to restrict to one spatial dimension. From now onwards, we consider D andD to be bounded Lipschitz domains in R and the spin accumulation s :D T → R 3 satisfies (2.3) and (2.4).

Regularity of weak solution
We wish to deduce improved stability properties for a weak solution of the problem (1.
i.e., there exists a constant C > 0 such that Recalling that H eff (m) = m and that we have assumed γ 0 = c = 1 in (1.1), we rewrite the (1.1) in the following form: where α 1 = α/(1 + α 2 ) and α 2 = 1/(1 + α 2 ), noting that |∇m| 2 = −m · m due to |m| = 1. This equivalence can be shown by taking the cross product to the left of both sides of (1.1), using |m| = 1 and the elementary identity a × (b × c) = b a, c − c a, b , and rearranging the resulting equation. Formally we multiply (4.1) with − m, and use the Cauchy-Schwarz inequality, the boundedness of m together with the Gagliardo-Nirenberg inequality for some σ,σ > 0 which can be chosen such that σ +σ < α 1 . Combining this with the estimate for s in (2.2) we conclude that m ∈ L 2 (H 2 ). Hence assertion i) follows. Proof of ii): First we note that (2.3) and (2.4) imply, after rearranging the equation, Then we formally multiply (4.2) by − s and integrate w.r.t x to obtain where By using the Cauchy-Schwarz inequality, the embedding H 1 → L ∞ , and (2.1), we have . Combining all the estimates we obtained from (4.3) By choosing sufficiently small and using the given assumption D * > 2βD * , we obtain after integrating over [0, t] and invoking Gronwall's Lemma In view of the estimation of i) and the estimate for s in (2.2), we complete the proof.
As we mentioned, we need improved regularity for the solution (m, s). To get improved regularity, we need to consider more regular control j ∈ L 2 (H 2 (D)).

)). Then the weak solution (m, s) of (1.1), (2.3) and (2.4) satisfies the following improved regularity:
Proof Proof of i): Let m 0 ∈ H 2 . Denote M = ∇m. Differentiating (4.1) with respect to the spatial variable, we formally have Multiplying (4.4) by − M (formally), using integration by parts, the Cauchy-Schwarz inequality, and the boundedness of m, we have after rearranging the equation for some small σ 1 > 0 and C ≡ C(σ 1 , α 1 , α 2 ) > 0 being a generic constant, so that by choosing σ 1 < α 1 and invoking Lemma 4.1, we deduce d dt By invoking Lemma 4.1 again and integrating over (0, t) we infer Gronwall's Lemma and Lemma 4.1 yield Again, we formally multiply (4.4) by ∂ t M, and use a similar argument along with the estimate (4.5) to conclude that Hence from (4.5) and (4.6) together with Lemma 4.
Proof of ii): Denote S = ∇s. Upon differentiating (2.3) with respect to the spatial variable, we formally see that S satisfies the following PDE: Multiplying (4.7) by − S (formally), and using the Cauchy-Schwarz inequality, the boundedness of m, and (2.1), we obtain By using the result in part i), we deduce (noting the assumption on D 0 ) We take = D * /2 and use Gronwall's lemma along with the assumptions on j to conclude Moreover, multiplying the equation (2.3) by ∂ t s, and then using the Cauchy-Schwarz inequality, part i) as well as (4.8), we have for any > 0 Since 0 < β < 1, we see that . Lastly, by formally multiplying (4.7) with ∂ t S, then applying the Cauchy-Schwarz inequality along with the estimate (4.8), one can arrive at the following estimate . This completes the proof.

Remark 4.1
In the proof of Lemma 4.2, we have used integration by parts formula, which may be made rigorous by a standard argument using the Faedo-Galerkin method which uses the related eigenfunctions of the given operator.

Optimization problem and its analysis
With the help of the Lagrange multiplier theorem ([18, Chapter 9, Theorem 1]), we now deduce the optimality system for the optimal solution of (1.6) where the con-straints (1.2) and (1.3) are replaced by (2.3) and (2.4). To this end, we define the spaces (4.9) Note that, in view of [10, Theorem 4, section 5.9.2, p. 306] we have The cost functional F defined in (1.6) could be re-interpreted as a function from M × S × J → R. With this set up, we now state the optimal control problem (1.6) in the following form.  X. An et al.

Proof
We calculate the directional derivatives of e 1 in the directions δm, δs and δj to obtain By using (4.10) we deduce Thus, e 1 is Gâteaux differentiable. Moreover, since e 1 is continuous at (m, s, j), the function e 1 is continuously Fréchet differentiable and the Fréchet derivative is given by (4.12). For the function e 2 , we have Moreover, we have the following bounds: Thus, e 2 is Gâteaux differentiable. Moreover, since e 2 is continuous at (m, s, j), the function e 2 is continuously Fréchet differentiable and the Fréchet derivative is given by (4.13). Similarly, we can show that e 3 and e 4 are continuously Fréchet differentiable, and their Fréchet derivatives are given by e 3 (m, s, j), (δm, δs, δj) = δm(0) and e 4 (m, s, j), (δm, δs, δj) = δs(0) .
Thus the function is continuously Fréchet differentiable. This completes the proof.
To apply the Lagrange multiplier theorem, one needs to check that a minimizer is a regular point of defined in (4.11). We recall that a point (m * , s * , j * ) is said to be a regular point of if e 1 (m * , s * , j * ), e 2 (m * , s * , j * ), e 3 (m * , s * , j * ), e 4 (m * , s * , j * ) are linearly independent. We have the following lemma. hold. We use several steps to solve the above coupled equations. Step (4.17) For f 1 and f 2 , this is not the case. Therefore, we define a projection operator in time.
To this end, let X be a Hilbert space with inner product (·, ·) X . For each time step k, we define the set For any f ∈ L 2 (0, T ; X), we define the projection k f ∈ P k via Existence of such projection k follows from the Lax-Milgram lemma. In view of the Cauchy-Schwarz inequality and (4.18) we see that and therefore we obtain for all f ∈ L 2 (0, T ; X). For any p k ∈ P k , we take φ k = k f − f + f − p k ∈ P k in (4.18) and use the Cauchy-Schwarz inequality to get Since p k ∈ P k is arbitrary and ∪ k>0 P k is dense in L 2 (0, T ; X), we get from above that (4.20) Let 1,k , 2,k , and 3,k denote the projections associated with X = L 2 , X = L 2 (D), and X = H 2 (D) respectively. In the following, we denote By using (4.19), we have Step II (Semi-discrete scheme and its solvability): We consider the following semi-discrete scheme for problems (4.14) and (4.15). Starting with m 0 = f 3 and Existence of m i +1 in step (i) given the existence of m i ∈ H 1 : Define a bilinear form A : One can use (4.17) for k is sufficiently small which can be chosen to be independent of the iteration step.
Here, in the last step we used again (4.17). Let Then it is linear and bounded, since Thus, by Lax-Milgram lemma, there exists a unique m i+1 ∈ H 1 such that (4.22) holds.

Existence of s i +1 in step (ii) given the existence of m i
.
We omit the details.
In order to derive the bound for m i from (4.26) and (4.28), we need to estimate . Choosing ψ = s i+1 as test function in (4.23) and using the boundedness of m * , the Cauchy-Schwarz inequality, and Young's inequality along with (2.1), we have for any θ > 0 By using (4.17), we obtain after rearranging the equation . By choosing θ such that D * − βD * − θ > 0, using (4.26) and (4.28), and multiplying the resulting equation by 2k, we deduce for any j > i. By summing over i from 0 to j − 1 for j ∈ {1, . . . , N}, we deduce By using (4.21), we have . By using the discrete Gronwall lemma, we deduce Using (4.29) in (4.26) and (4.28), we obtain max 0≤i≤N where C is a constant independent of the time step size k. Hence, by choosing ψ = − s i+1 in (4.23) and using integration by parts along with (4.17), we have Note, in the calculation above, we used the boundedness of D 0 (D 0 ∈ H 2 (D)) to have ∇D 0 L ∞ ≤ C. (4.31) By multiplying the above equation by 2k, using (4.29) and (4.30), we deduce, after rearranging the equation and choosing such that D * − − βD * > 0, By summing over i from 0 to j − 1 for j ∈ {1, . . . , N}, using (4.21), and invoking the discrete Gronwall lemma, we deduce for some constant C > 0 independent of k. Combining this equation with (4.29) gives Step IV (Continuation and its bound): Moreover, we define the piecewise constant interpolants (in time) X + k and X − k as follows: Choosing the test function φ = d t m i+1 in (4.22) and using integration by parts, we obtain (recalling (4.17) again) Hence, by using (4.28) and (4.32), as well as rearranging the equation, we have Integrating over the time interval (t i , t i+1 ) and summing over i, we obtain where we used (4.21) and (4.28) in the last step. Next we show the boundedness of the sequence . Then using the test function ψ = d t s i+1 in (4.23) and integration by parts, we get (recalling (4.17) and (4.31)) or equivalently, Integrating over the time interval (t i , t i+1 ) and summing over i, we get where in the last step we used (4.21) and (4.32).
Step V (Solvability of (4.14) and (4.15)): In view of (4.30), (4.32), (4.34), and (4.35), there exists a constant C > 0, independent of the discretization parameter k > 0, such that Inequality (4.30) and simple calculations reveal D)). Therefore, it follows from (4.36) that there exists (m, s) ∈ M × S (defined by (4.9)) such that the following statements hold: (H 1 (D)), S ± k → s in L 2 (L 2 (D)). (4.37) Moreover, if M * ± k are defined from m * i = m * (t i ) by (4.33) then, since m * ∈ C([0, T ], H 2 ), it can be shown, using the uniform continuity of m * , that We now prove that (m, s) given in (4.37) is a solution to (4.14) and (4.15). It follows from (4.22) and (4.23) that M k , M k , M ± k , S k , and S ± k satisfy the following equations We discuss the passing to the limit (when k → 0) of (4.39) only. Similar arguments hold for (4.40). The convergence in (4.37) implies, for any φ ∈ C ∞ (C ∞ (D)), The convergence of the linear term can be obtained by using (4.20) as follows: The convergence of the nonlinear terms in (4.39) can be derived by using (4.36), (4.37), and (4.38). We present detailed estimates for one term, namely The convergence of other nonlinear terms in (4.39) can be shown in the same manner. Therefore, by letting k → 0 in (4.39) we deduce or equivalently (by using integration by parts), which implies the first equation in (4.14).
To show the second equation in (4.14), namely m(0) = f 3 , we choose the test function φ ∈ C ∞ (C ∞ (D)) such that φ(T ) = 0. Then integration by parts gives X. An et al.
Using integration by parts again gives the required result. Similarly, one can easily show that (m, s) satisfies (4.15), completing the proof of the lemma.

Regularity of the adjoint and optimal control
In view of the Lagrange multiplier theorem, the adjoint variables z 1 and z 2 satisfy (z 1 , z 2 ) ∈ L 2 (L 2 ) × L 2 (L 2 (D)). But one can expect better regularity properties for these adjoint variables and the optimal control j * . Regarding this, we have the following lemma.
Thus, we complete the proof of assertion i). Proof of ii): Test (2.5) with j * and − j * and then integrate with respect to x along with the Cauchy-Schwarz inequality. The resulting inequality becomes Thanks to i), we see from the above estimate that j * ∈ L ∞ (H 3 (D)). Proof of iii): Testing (2.6) with z 1 , and − z 1 , we have, thanks to the Cauchy-Schwarz inequality, the boundedness of m * , and the embedding of L ∞ between L 2 and H 1 for some > 0. We now use part i) and part ii) above, Lemma 4.2, and the condition on z 1 (T ) given in (2.8) to conclude that z 1 ∈ L ∞ (H 1 ) ∩ L 2 (H 2 ). As in the proof of part i), by testing (2.6) with ∂ t z 1 , one can easily show that ∂ t z 1 2 L 2 (L 2 ) ≤ C, which completes the proof of iii).

Numerical experiments
The following numerical experiments are carried out using the FEniCS automated code generation system (https://fenicsproject.org) version 2019 and the tlm adjoint library (https://github.com/jrmaddison/tlm adjoint). FEniCS enables a high-level syntax representation of the complex numerical equations, as well as an efficient implementation using automated tools. On the other hand, the tlm adjoint library is used to derive the associated adjoint model, which computes the required derivatives of the cost functional [11,15,19,27]. In general, the implementation only requires the users to define a suitable finite element function space, the computational domain and mesh, the weak formulation, and a few other specifics. The FEniCS system and the tlm adjoint library enable automatic derivation of the finite element equation, the discrete tangent-linear and adjoint models. Nevertheless, we would like to lay out more details for pedagogical purpose. Generally, the PDE constrained optimization algorithm contains two main parts. The first part is to solve the forward problem, while the second part is to achieve optimization.
We elaborate the first part here. First we note that the weak equations in Definition 2.1 can be rewritten as where we recover the constants γ 0 , c, β , e, μ β λ 1 , and λ 2 as in problem (1.1)-(1.3). The fully discrete schemes to solve these equations are (5.1) and (5.2) below. For the time discretization, we recall (4.16) and use the backward Euler scheme. The spatial discretization is determined by a shape-regular triangulation T h of D andT h ofD into tetrahedra such that the two triangulations agree on D.
The finite element spaces are defined by 3 for all τ ∈T h }, d (i) in which the functional F is decreasing. We used the L-BFGS-B algorithm [22], which is an iterative method for solving unconstrained nonlinear optimization problems from the Scipy library [32]. (e) Set i = i + 1 and return to (a).
Step (d) is crucial as the main task of the gradient based optimization is to determine an improved choice of control variable in order to achieve fast convergence to the optimal solution.
In Section 5.1 we show the numerical results for the case of one-dimensional domains D andD. Even though the analysis presented in previous sections does not apply for multi-dimensional domains, we still experiment with D =D ⊂ R 3 ; the results are shown in Section 5.2. In all experiments in both sections, the cost functional F (π), see (1.5), is defined with the terminal payoff
The stopping criterion is We apply uniform partitions with space mesh size h = 0.01 and time step k = 0.0005. We observe that when a larger value of k is chosen, e.g., k = 0.5, the numerical simulation becomes unstable and Newton's method used to solve the nonlinear system (5.1) fails to converge.

No targets
Example 1: In this example, we show the evolutions of magnetization m and the spin accumulation s when the control variable j = 0 and there are no targets (Fig. 1), with

Single targets
Example 2: In this example, the initial states, total running time, and the assimilation target are  The simulation is presented in Fig. 3 (1a)-(1e). The simulation is presented in Fig. 4 (1a)-(1e). For all three examples 2, 3, and 4, we show in Fig. 5 two comparisons. The left column (a) plots the cost functional F versus the optimization iteration. The right column (b) plots the relative change δF /F versus the optimization iteration. A general decreasing trend can be seen which indicates the numerical method is valid. (5.5) Figure 6 shows the evolution of the optimal m * with respect to time t. Figure 6 (1d) is the simulation at the switching moment t = 0.3.

Multiple cost functional settings
In  The assimilation process for Example 6 is shown in Fig. 7 and for Example 7 in Fig. 8.
The same stopping criterion (5.3) is used. We apply uniform partition with space mesh size h = 0. The progress of the assimilation is presented in Fig. 9 (a)-(d). It is noted that in order to achieve optimal control in this example, we need to add an extra external field, namely the Zeeman field, to the effective field in (1.4). Otherwise, the minimization of the cost functional cannot be obtained. The effective field in this example becomes where e 3 = (0, 0, 1) and c is some physical constant, chosen to be 1 in this example. The assimilation is shown in Fig. 10  While the time step is chosen to be k = 0.01 as in the previous examples, here we have to choose a smaller spatial step, namely h = 0.05. This is due to the nonconstant initial data m 0 which results in larger interpolation errors at non-nodal points in D, even though at the nodal points normalization has been performed to ensure that |m| = 1. These errors can be observed in Fig. 11 (c) (red color area). The progress of the assimilation is shown in Fig. 11 (a)-(d).

Discussion
The numerical simulations in Sections 5.1 and 5.2 showcase the availability of optimal control for a coupled spin-polarized current and magnetization system in a simple geometry. It is noted that, in order to provide a realistic physical interpretation and insight, a couple of issues need to be addressed. First, a more sophisticated geometry is needed. The typical case is a multilayer structure that consists of two ferromagnetic layers of given thickness that are separated by a nonmagnetic layer [1,2,28].
Second, more external fields should be included to the effective fields in the LLG equation. However, because this paper focuses on providing a mathematical framework for the optimal control problem, we only use the simple geometry to validate the mathematical analysis. We plan to investigate the optimal control of the coupled system in a more complicated setting in the near future.

Concluding remarks
In this paper, we proved the existence of the optimal solution of a coupled spin drift-diffusion Landau-Lifshitz-Gilbert system. We also showed the existence of the adjoint variables which define the first-order optimality system to be solved for the

Conflict of interest Authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.