Yuto Miyatake

Yuto Miyatake, David Cohen, Daisuke Furihata et al.

We present novel geometric numerical integrators for Hunter--Saxton-like equations by means of new multi-symplectic formulations and known Hamiltonian structures of the problems. We consider the Hunter--Saxton equation, the modified Hunter--Saxton equation, and the two-component Hunter--Saxton equation. Multi-symplectic discretisations based on these new formulations of the problems are exemplified by means of the explicit Euler box scheme, and Hamiltonian-preserving discretisations are exemplified by means of the discrete variational derivative method. We explain and justify the correct treatment of boundary conditions in a unified manner. This is necessary for a proper numerical implementation of these equations and was never explicitly clarified in the literature before, to the best of our knowledge. Finally, numerical experiments demonstrate the favourable behaviour of the proposed numerical integrators.

Yuto Miyatake

Since the expense of the numerical integration of large scale dynamical systems is often computationally prohibitive, model reduction methods, which approximate such systems by simpler and much lower order ones, are often employed to reduce the computational effort. In this paper, for dynamical systems with a first integral, new structure-preserving model reduction approaches are presented that yield reduced-order systems while preserving the first integral. We apply energy-preserving integrators to the reduced-order systems and show some numerical experiments that demonstrate the favourable behaviour of the proposed approaches.

Yuto Miyatake, Tomohiro Sogabe, Shao-Liang Zhang

Because the expense of estimating the optimal value of the relaxation parameter in the successive over-relaxation (SOR) method is usually prohibitive, the parameter is often adaptively controlled. In this paper, new adaptive SOR methods are presented that are applicable to a variety of symmetric positive definite linear systems and do not require additional matrix-vector products when updating the parameter. To this end, we regard the SOR method as an algorithm for minimising a certain objective function, which yields an interpretation of the relaxation parameter as the step size following a certain change of variables. This interpretation enables us to adaptively control the step size based on some line search techniques, such as the Wolfe conditions. Numerical examples demonstrate the favourable behaviour of the proposed methods.

Masaya Oozawa, Tomohiro Sogabe, Yuto Miyatake et al.

We consider the T-congruence Sylvester equation $AX+X^{\rm T}B=C$, where $A\in \mathbb R^{m\times n}$, $B\in \mathbb R^{n\times m}$ and $C\in \mathbb R^{m\times m}$ are given, and matrix $X \in \mathbb R^{n\times m}$ is to be determined. The T-congruence Sylvester equation has recently attracted attention because of a relationship with palindromic eigenvalue problems. For example, necessary and sufficient conditions for the existence and uniqueness of solutions, and numerical solvers have been intensively studied. In this note, we will show that, under a certain condition, the T-congruence Sylvester equation can be transformed into the Lyapunov equation. This may lead to further properties and efficient numerical solvers by utilizing a great deal of studies on the Lyapunov equation.

Fuminori Tatsuoka, Tomohiro Sogabe, Yuto Miyatake et al.

Incremental Newton (IN) iteration, proposed by Iannazzo, is stable for computing the matrix $p$th root, and its computational cost is $\mathcal{O}(n^3p)$ flops per iteration. In this paper, a cost-efficient variant of IN iteration is presented. The computational cost of the variant well agrees with $\mathcal{O} (n^3 \log p)$ flops per iteration, if $p$ is up to at least 100.

Yuto Miyatake, Geonsik Eom, Tomohiro Sogabe et al.

We consider the numerical integration of the Hunter--Saxton equation, which models the propagation of weakly nonlinear orientation waves. For the equation, we present two weak forms and their Galerkin discretizations. The Galerkin schemes preserve the Hamiltonian of the equation and can be implemented with cheap $H^1$ elements. Numerical experiments confirm the effectiveness of the schemes.

Fuminori Tatsuoka, Yuto Miyatake, Tomohiro Sogabe

Several methods for computing the action of the matrix exponential $\mathrm{e}^{\boldsymbol{A}} \boldsymbol{b}$ are expressed by substituting $\boldsymbol{A}$ into a rational approximation of the scalar exponential function. The error of such methods can be estimated using the numerical range of $\boldsymbol{A}$, which enables the computation of $\mathrm{e}^{\boldsymbol{A}}\boldsymbol{b}$ with a prescribed accuracy. However, when the input matrix has the structure $\boldsymbol{A} = τ\boldsymbol{M}^{-1} \boldsymbol{K}$, this approach is challenging because computing the bounding box of numerical range is difficult and the numerical range may be too large to construct rational approximations on it. In this paper, focusing on the case where $\boldsymbol{M}$ is a well-conditioned symmetric positive definite matrix, we propose considering the numerical range of a similarity transformed matrix of $\boldsymbol{A}$. The numerical range of transformed matrix is not only numerically computable but can also be theoretically bounded depending on properties of $\boldsymbol{K}$. Numerical experiments confirm that the computations can be performed within the prescribed error tolerance.

Hiroki Takeda, Yuto Miyatake, Daisuke Furihata

Tensor Train (TT) decomposition is a powerful technique for analyzing high-dimensional data. Existing algorithms for computing TT decompositions can be categorized into two main types: conventional batch-based approaches and recursive online methods. In the context of streaming data, batch methods typically achieve higher reconstruction accuracy but often suffer from memory exhaustion, while online methods provide greater computational efficiency. In this work, we introduce Online TT-ALS (Alternating Least Squares), an algorithm that sequentially enforces orthogonality constraints. This approach allows for efficient and exact updates of the core tensor while maintaining high reconstruction accuracy. Theoretically, we prove that enforcing these orthogonal gauge constraints guarantees monotonic decrease of the local objective function and temporal smoothness. Computationally, our deterministic single-sweep update reduces the rank dependence from quadratic to linear, achieving an overall complexity of $\mathcal{O}(I^{n-1} r)$. Experimental results demonstrate that the proposed method outperforms existing online techniques not only in terms of mathematical approximation accuracy but also in human perception-based video quality metrics. Furthermore, compared to recent deep learning-based paradigms, our algebraic approach achieves speedups of several orders of magnitude. Consequently, our method exhibits high computational efficiency and is suitable for low-latency real-time processing applications.

Yuto Miyatake

In the study of energy-preserving methods for Hamiltonian systems, polynomial continuous-stage Runge--Kutta methods play an important role. Necessary and sufficient conditions for such methods to be energy-preserving have already been established. They are energy-preserving if the matrix $M\in \mathbb{R}^{s\times s}$ defining the method is symmetric, and the converse holds under the assumption that a certain $s\times \infty$ matrix $Φ^\mathrm{CSRK}$ has full row rank. It was conjectured in Remark 3 in Miyatake and Butcher (SIAM J. Numer. Anal., 2016) that the full-rank assumption should always hold for every consistent polynomial continuous-stage Runge--Kutta method. In this paper, we prove the conjecture by showing that the matrix $Φ^\mathrm{CSRK}$ has full row rank under the standard consistency condition. The proof is a direct application of the polynomial moment problem solved by Pakovich and Muzychuk (Proc. Lond. Math. Soc., 2009).

Haoyu Chu, Yuto Miyatake, Wenjun Cui et al.

Recently, there has been growing interest in using physics-informed neural networks (PINNs) to solve differential equations. However, the preservation of structure, such as energy and stability, in a suitable manner has yet to be established. This limitation could be a potential reason why the learning process for PINNs is not always efficient and the numerical results may suggest nonphysical behavior. Besides, there is little research on their applications on downstream tasks. To address these issues, we propose structure-preserving PINNs to improve their performance and broaden their applications for downstream tasks. Firstly, by leveraging prior knowledge about the physical system, a structure-preserving loss function is designed to assist the PINN in learning the underlying structure. Secondly, a framework that utilizes structure-preserving PINN for robust image recognition is proposed. Here, preserving the Lyapunov structure of the underlying system ensures the stability of the system. Experimental results demonstrate that the proposed method improves the numerical accuracy of PINNs for partial differential equations. Furthermore, the robustness of the model against adversarial perturbations in image data is enhanced.

Kosuke Akita, Yuto Miyatake, Daisuke Furihata

In data assimilation, state estimation is not straightforward when the observation operator is unknown. This study proposes a method for composing a surrogate operator when the true operator is unknown. A neural network is used to improve the surrogate model iteratively to decrease the difference between the observations and the results of the surrogate model. A twin experiment suggests that the proposed method outperforms approaches that tentatively use a specific operator throughout the data assimilation process.

Takashi Matsubara, Yuto Miyatake, Takaharu Yaguchi

A neural network model of a differential equation, namely neural ODE, has enabled the learning of continuous-time dynamical systems and probabilistic distributions with high accuracy. The neural ODE uses the same network repeatedly during a numerical integration. The memory consumption of the backpropagation algorithm is proportional to the number of uses times the network size. This is true even if a checkpointing scheme divides the computation graph into sub-graphs. Otherwise, the adjoint method obtains a gradient by a numerical integration backward in time. Although this method consumes memory only for a single network use, it requires high computational cost to suppress numerical errors. This study proposes the symplectic adjoint method, which is an adjoint method solved by a symplectic integrator. The symplectic adjoint method obtains the exact gradient (up to rounding error) with memory proportional to the number of uses plus the network size. The experimental results demonstrate that the symplectic adjoint method consumes much less memory than the naive backpropagation algorithm and checkpointing schemes, performs faster than the adjoint method, and is more robust to rounding errors.

Kimiaki Kinugasa, Yuto Miyatake, Takayasu Matsuo

The numerical integration of the Benjamin and Benjamin--Ono equations are considered. They are non-local partial differential equations involving the Hilbert transform, and due to this, so far quite few structure-preserving integrators have been proposed. In this paper, a new reformulation of the equations is stated, and new structure-preserving discretizations are proposed based on it. Numerical experiments confirm the effectiveness of the proposed integrators.