Applied Mathematics and Computational Science Program
For more information visit: https://cemse.kaust.edu.sa/amcs
Recent Submissions

Impact of Evanescence Process on ThreeDimensional SubDiffusion based Molecular Communication Channel(IEEE Transactions on NanoBioscience, Institute of Electrical and Electronics Engineers (IEEE), 20230321) [Article]In most of the existing works of molecular communication (MC), the standard diffusion environment is taken into account where the mean square displacement (MSD) of an information molecule (IM) scales linearly with time. On the contrary, this work considers the subdiffusion motion that appears in crowded and complex (porous or fractal) environments (movement of the particles in the living cells) where the particle’s MSD scales as a fractional order power law in time. Moreover, we examine an additional evanescence process resulting from which the molecules can degrade before hitting the boundary of the receiver (RX). Thus, in this work, we present a 3D MC system with a point transmitter (TX) and the spherical RX with the subdiffusive behavior of an IM along with its evanescence. Furthermore, an IM’s closedform expressions for the arrival probability and the first passage time density (FPTD) are emulated in the above context. Additionally, we investigate the performance of MC by using the concentrationbased modulation technique in a subdiffusion channel. Finally, the considered MC channel is exploited in terms of the probability of detection, probability of false alarm, and probability of error for different parameters such as the reaction rate, fractional power, and radius of the RX.

Finite element discretization of a biological network formation system: a preliminary study(arXiv, 20230319) [Preprint]A finite element discretization is developed for the CaiHu model, describing the formation of biological networks. The model consists of a non linear elliptic equation for the pressure p and a non linear reactiondiffusion equation for the conductivity tensor C. The problem requires high resolution due to the presence of multiple scales, the stiffness in all its components and the non linearities. We propose a low order finite element discretization in space coupled with a semiimplicit time advancing scheme. The code is validated with several numerical tests performed with various choices for the parameters involved in the system. In absence of the exact solution, we apply Richardson extrapolation technique to estimate the order of the method.

Generalization of the Orthodiagonal Involutive Type of Kokotsakis Flexible Polyhedra(arXiv, 20230319) [Preprint]In this paper we introduce and study a remarkable class of mechanisms formed by a 3×3 arrangement of rigid and skew quadrilateral faces with revolute joints at the common edges. These Kokotsakistype mechanisms with a quadrangular base and nonplanar faces are a generalization of Izmestiev's orthodiagonal involutive type of Kokotsakis polyhedra formed by planar quadrilateral faces. Our algebraic approach yields a complete characterization of all complexes of the orthodiagonal involutive type. It is shown that one has 8 degrees of freedom to construct such mechanisms. This is illustrated by several examples, including cases that are not possible with planar faces.

Test and Visualization of Covariance Properties for Multivariate SpatioTemporal Random Fields(Journal of Computational and Graphical Statistics, Informa UK Limited, 20230316) [Article]The prevalence of multivariate spacetime data collected from monitoring networks and satellites, or generated from numerical models, has brought much attention to multivariate spatiotemporal statistical models, where the covariance function plays a key role in modeling, inference, and prediction. For multivariate spacetime data, understanding the spatiotemporal variability, within and across variables, is essential in employing a realistic covariance model. Meanwhile, the complexity of generic covariances often makes model fitting very challenging, and simplified covariance structures, including symmetry and separability, can reduce the model complexity and facilitate the inference procedure. However, a careful examination of these properties is needed in real applications. In the work presented here, we formally define these properties for multivariate spatiotemporal random fields and use functional data analysis techniques to visualize them, hence providing intuitive interpretations. We then propose a rigorous rankbased testing procedure to conclude whether the simplified properties of covariance are suitable for the underlying multivariate spacetime data. The good performance of our method is illustrated through synthetic data, for which we know the true structure. We also investigate the covariance of bivariate wind speed, a key variable in renewable energy, over a coastal and an inland area in Saudi Arabia. The Supplementary Material is available online, including the R code for our developed methods.

Semantic Segmentation of Mesoscale Eddies in the Arabian Sea: A Deep Learning Approach(Remote Sensing, MDPI AG, 20230310) [Article]Detecting mesoscale ocean eddies provides a better understanding of the oceanic processes that govern the transport of salt, heat, and carbon. Established eddy detection techniques rely on physical or geometric criteria, and they notoriously fail to predict eddies that are neither circular nor elliptical in shape. Recently, deep learning techniques have been applied for semantic segmentation of mesoscale eddies, relying on the outputs of traditional eddy detection algorithms to supervise the training of the neural network. However, this approach limits the network’s predictions because the available annotations are either circular or elliptical. Moreover, current approaches depend on the seasurface height, temperature, or currents as inputs to the network, and these data may not provide all the information necessary to accurately segment eddies. In the present work, we have trained a neural network for the semantic segmentation of eddies using humanbased—and expertvalidated—annotations of eddies in the Arabian Sea. Training with humanannotated datasets enables the network predictions to include more complex geometries, which occur commonly in the real ocean. We then examine the impact of different combinations of input surface variables on the segmentation performance of the network. The results indicate that providing additional surface variables as inputs to the network improves the accuracy of the predictions by approximately 5%. We have further finetuned another pretrained neural network to segment eddies and achieved a reduced overall training time and higher accuracy compared to the results from a network trained from scratch.

Modeling metallic fatigue data using the BirnbaumSaunders distribution(arXiv, 20230309) [Preprint]This work employs the BirnbaumSaunders distribution to model metallic fatigue and compares its performance to fatiguelimit models based on the normal distribution. First, we fit data for 85 fatigue experiments with constant amplitude cyclic loading applied to unnotched sheet specimens of 75ST6 aluminum alloys. The fit obtained by the BirnbaumSaunders distribution is noticeably better than the normal distribution. Then, we define new equivalent stress for two fatigue experiment types: tensioncompression and tensiontension. With the new equivalent stress, the statistical fit improves for both distributions, with a slight preference for the BirnbaumSaunders distribution. In addition, we analyze a dataset of rotatingbending fatigue experiments applied to 101 round bar specimens of 75ST6 aluminum. Finally, we consider a wellknown dataset of bending tests applied to 125 specimens of carbon laminate. Overall, the BirnbaumSaunders distribution provides better fit results under fatiguelimit models with various experimental setups.

Analysis of a Class of Multilevel Markov Chain Monte Carlo Algorithms Based on Independent Metropolis–Hastings(SIAM/ASA Journal on Uncertainty Quantification, Society for Industrial & Applied Mathematics (SIAM), 20230303) [Article]In this work, we present, analyze, and implement a class of multilevel Markov chain Monte Carlo (MLMCMC) algorithms based on independent Metropolis–Hastings proposals for Bayesian inverse problems. In this context, the likelihood function involves solving a complex differential model, which is then approximated on a sequence of increasingly accurate discretizations. The key point of this algorithm is to construct highly coupled Markov chains together with the standard multilevel Monte Carlo argument to obtain a better costtolerance complexity than a singlelevel MCMC algorithm. Our method extends the ideas of Dodwell et al., [SIAM/ASA J. Uncertain. Quantif., 3 (2015), pp. 1075–1108] to a wider range of proposal distributions. We present a thorough convergence analysis of the MLMCMC method proposed, and show, in particular, that (i) under some mild conditions on the (independent) proposals and the family of posteriors, there exists a unique invariant probability measure for the coupled chains generated by our method, and (ii) that such coupled chains are uniformly ergodic. We also generalize the costtolerance theorem of Dodwell et al. to our wider class of MLMCMC algorithms. Finally, we propose a selftuning continuationtype MLMCMC algorithm. The presented method is tested on an array of academic examples, where some of our theoretical results are numerically verified. These numerical experiments evidence how our extended MLMCMC method is robust when targeting some pathological posteriors, for which some of the previously proposed MLMCMC algorithms fail.

Doubleloop quasiMonte Carlo estimator for nested integration(arXiv, 20230227) [Preprint]Nested integration arises when a nonlinear function is applied to an integrand, and the result is integrated again, which is common in engineering problems, such as optimal experimental design, where typically neither integral has a closedform expression. Using the Monte Carlo method to approximate both integrals leads to a doubleloop Monte Carlo estimator, which is often prohibitively expensive, as the estimation of the outer integral has bias relative to the variance of the inner integrand. For the case where the inner integrand is only approximately given, additional bias is added to the estimation of the outer integral. Variance reduction methods, such as importance sampling, have been used successfully to make computations more affordable. Furthermore, random samples can be replaced with deterministic lowdiscrepancy sequences, leading to quasiMonte Carlo techniques. Randomizing the lowdiscrepancy sequences simplifies the error analysis of the proposed doubleloop quasiMonte Carlo estimator. To our knowledge, no comprehensive error analysis exists yet for truly nested randomized quasiMonte Carlo estimation (i.e., for estimators with lowdiscrepancy sequences for both estimations). We derive asymptotic error bounds and a method to obtain the optimal number of samples for both integral approximations. Then, we demonstrate the computational savings of this approach compared to standard nested (i.e., doubleloop) Monte Carlo integration when estimating the expected information gain via two examples from Bayesian optimal experimental design, the latter of which involves an experiment from solid mechanics.

A mode of convergence arising in diffusive relaxation(arXiv, 20230227) [Preprint]In this work, a mode of convergence for measurable functions is introduced. A related notion of Cauchy sequence is given and it is proved that this notion of convergence is complete in the sense that Cauchy sequences converge. Moreover, the preservation of convergence under composition is investigated. The origin of this mode of convergence lies in the path of proving that the density of a Euler system converges almost everywhere (up to subsequences) towards the density of a nonlinear diffusion system, as a consequence of the convergence in the relaxation limit.

Existence and weakstrong uniqueness for MaxwellStefanCahnHilliard systems(Accepted by Annales de l'Institut Henri Poincaré C, Analyse non linéaire, Elsevier, 20230226) [Article]A Maxwell–Stefan system for fluid mixtures with driving forces depending on Cahn–Hilliardtype chemical potentials is analyzed. The corresponding parabolic crossdiffusion equations contain fourthorder derivatives and are considered in a bounded domain with noflux boundary conditions. The nonconvex part of the energy is assumed to have a bounded Hessian. The main difficulty of the analysis is the degeneracy of the diffusion matrix, which is overcome by proving the positive definiteness of the matrix on a subspace and using the Bott–Duffin matrix inverse. The global existence of weak solutions and a weakstrong uniqueness property are shown by a careful combination of (relative) energy and entropy estimates, yielding H2 (Ω) bounds for the densities, which cannot be obtained from the energy or entropy inequalities alone.

Datadriven uncertainty quantification for constrained stochastic differential equations and application to solar photovoltaic power forecast data(arXiv, 20230225) [Preprint]In this work, we extend the datadriven It\^{o} stochastic differential equation (SDE) framework for the pathwise assessment of shortterm forecast errors to account for the timedependent upper bound that naturally constrains the observable historical data and forecast. We propose a new nonlinear and timeinhomogeneous SDE model with a Jacobitype diffusion term for the phenomenon of interest, simultaneously driven by the forecast and the constraining upper bound. We rigorously demonstrate the existence and uniqueness of a strong solution to the SDE model by imposing a condition for the timevarying meanreversion parameter appearing in the drift term. The normalized forecast function is thresholded to keep such meanreversion parameters bounded. The SDE model parameter calibration also covers the thresholding parameter of the normalized forecast by applying a novel iterative twostage optimization procedure to userselected approximations of the likelihood function. Another novel contribution is estimating the transition density of the forecast error process, not known analytically in a closed form, through a tailored kernel smoothing technique with the control variate method. We fit the model to the 2019 photovoltaic (PV) solar power daily production and forecast data in Uruguay, computing the daily maximum solar PV production estimation. Two statistical versions of the constrained SDE model are fit, with the beta and truncated normal distributions as proxies for the transition density. Empirical results include simulations of the normalized solar PV power production and pathwise confidence bands generated through an indirect inference method. An objective comparison of optimal parametric points associated with the two selected statistical approximations is provided by applying the innovative kernel density estimation technique of the transition function of the forecast error process.

Reshaping Geostatistical Modeling and Prediction for ExtremeScale Environmental Applications(IEEE, 20230223) [Conference Paper]We extend the capability of spacetime geostatistical modeling using algebraic approximations, illustrating applicationexpected accuracy worthy of double precision from majority lowprecision computations and lowrank matrix approximations. We exploit the mathematical structure of the dense covariance matrix whose inverse action and determinant are repeatedly required in Gaussian loglikelihood optimization. Geostatistics augments firstprinciples modeling approaches for the prediction of environmental phenomena given the availability of measurements at a large number of locations; however, traditional Choleskybased approaches grow cubically in complexity, gating practical extension to continental and global datasets now available. We combine the linear algebraic contributions of mixedprecision and lowrank computations within a tile based Cholesky solver with ondemand casting of precisions and dynamic runtime support from PaRSEC to orchestrate tasks and data movement. Our adaptive approach scales on various systems and leverages the Fujitsu A64FX nodes of Fugaku to achieve up to 12X performance speedup against the highly optimized dense Cholesky implementation.

Physicsinformed Spectral Learning: the Discrete HelmholtzHodge Decomposition(arXiv, 20230221) [Preprint]In this work, we further develop the Physicsinformed Spectral Learning (PiSL) by Espath et al. \cite{Esp21} based on a discrete L2 projection to solve the discrete HodgeHelmholtz decomposition from sparse data. Within this physicsinformed statistical learning framework, we adaptively build a sparse set of Fourier basis functions with corresponding coefficients by solving a sequence of minimization problems where the set of basis functions is augmented greedily at each optimization problem. Moreover, our PiSL computational framework enjoys spectral (exponential) convergence. We regularize the minimization problems with the seminorm of the fractional Sobolev space in a Tikhonov fashion. In the Fourier setting, the divergence and curlfree constraints become a finite set of linear algebraic equations. The proposed computational framework combines supervised and unsupervised learning techniques in that we use data concomitantly with the projection onto divergence and curlfree spaces. We assess the capabilities of our method in various numerical examples including the `Storm of the Century' with satellite data from 1993.

An Improved Unbiased Particle Filter(arXiv, 20230220) [Preprint]In this paper we consider the filtering of partially observed multidimensional diffusion processes that are observed regularly at discrete times. We assume that, for numerical reasons, one has to timediscretize the diffusion process which typically leads to filtering that is subject to discretization bias. The approach in [16] establishes that when only having access to the timediscretized diffusion it is possible to remove the discretization bias with an estimator of finite variance. We improve on the method in [16] by introducing a modified estimator based on the recent work of [17]. We show that this new estimator is unbiased and has finite variance. Moreover, we conjecture and verify in numerical simulations that substantial gains are obtained. That is, for a given mean square error (MSE) and a particular class of multidimensional diffusion, the cost to achieve the said MSE falls.

DISPERSIVE SHOCKS IN DIFFUSIVEDISPERSIVE APPROXIMATIONS OF ELASTICITY AND QUANTUMHYDRODYNAMICS(QUARTERLY OF APPLIED MATHEMATICS, American Mathematical Society (AMS), 20230217) [Article]The aim is to assess the combined effect of diffusion and dispersion on shocks in the moderate dispersion regime. For a diffusive dispersive approximation of the equations of onedimensional elasticity (or psystem), we study convergence of traveling waves to shocks. The problem is recast as a Hamiltonian system with small friction, and an analysis of the length of oscillations yields convergence in the moderate dispersion regime ε, δ → 0 with δ = o(ε), under hypotheses that the limiting shock is admissible according to the Liu Econdition and is not a contact discontinuity at either end state. A similar convergence result is proved for traveling waves of the quantum hydrodynamic system with artificial viscosity as well as for a viscous PeregrineBoussinesq system where traveling waves model undular bores, in all cases in the moderate dispersion regime.

Bounds on Depth of Decision Trees Derived from Decision Rule Systems(arXiv, 20230214) [Preprint]Systems of decision rules and decision trees are widely used as a means for knowledge representation, as classifiers, and as algorithms. They are among the most interpretable models for classifying and representing knowledge. The study of relationships between these two models is an important task of computer science. It is easy to transform a decision tree into a decision rule system. The inverse transformation is a more difficult task. In this paper, we study unimprovable upper and lower bounds on the minimum depth of decision trees derived from decision rule systems depending on the various parameters of these systems.

PhysicsInformed Deep Neural Network for BackwardinTime Prediction: Application to Rayleigh–Bénard Convection(Artificial Intelligence for the Earth Systems, American Meteorological Society, 20230214) [Article]Backwardintime predictions are needed to better understand the underlying dynamics of physical fluid flows and improve future forecasts. However, integrating fluid flows backward in time is challenging because of numerical instabilities caused by the diffusive nature of the fluid systems and nonlinearities of the governing equations. Although this problem has been long addressed using a nonpositive diffusion coefficient when integrating backward, it is notoriously inaccurate. In this study, a physicsinformed deep neural network (PIDNN) is presented to predict past states of a dissipative dynamical system from snapshots of the subsequent evolution of the system state. The performance of the PIDNN is investigated using several systematic numerical experiments and the accuracy of the backwardintime predictions is evaluated in terms of different error metrics. The proposed PIDNN can predict the previous state of the Rayleigh–Bénard convection with an 8time step average normalized ℓ2error of less than 2% for a turbulent flow at a Rayleigh number of 105.

Decision Trees for Binary SubwordClosed Languages(Entropy, MDPI AG, 20230214) [Article]In this paper, we study arbitrary subwordclosed languages over the alphabet {0,1} (binary subwordclosed languages). For the set of words L(n) of the length n belonging to a binary subwordclosed language L, we investigate the depth of the decision trees solving the recognition and the membership problems deterministically and nondeterministically. In the case of the recognition problem, for a given word from L(n), we should recognize it using queries, each of which, for some i∈{1,…,n}, returns the ith letter of the word. In the case of the membership problem, for a given word over the alphabet {0,1} of the length n, we should recognize if it belongs to the set L(n) using the same queries. With the growth of n, the minimum depth of the decision trees solving the problem of recognition deterministically is either bounded from above by a constant or grows as a logarithm, or linearly. For other types of trees and problems (decision trees solving the problem of recognition nondeterministically and decision trees solving the membership problem deterministically and nondeterministically), with the growth of n, the minimum depth of the decision trees is either bounded from above by a constant or grows linearly. We study the joint behavior of the minimum depths of the considered four types of decision trees and describe five complexity classes of binary subwordclosed languages.

Uncertainty quantification in coastal aquifers using the multilevel Monte Carlo method(arXiv, 20230213) [Preprint]We consider a class of densitydriven flow problems. We are particularly interested in the problem of the salinization of coastal aquifers. We consider the Henry saltwater intrusion problem with uncertain porosity, permeability, and recharge parameters as a test case. The reason for the presence of uncertainties is the lack of knowledge, inaccurate measurements, and inability to measure parameters at each spatial or time location. This problem is nonlinear and timedependent. The solution is the salt mass fraction, which is uncertain and changes in time. Uncertainties in porosity, permeability, recharge, and mass fraction are modeled using random fields. This work investigates the applicability of the wellknown multilevel Monte Carlo (MLMC) method for such problems. The MLMC method can reduce the total computational and storage costs. Moreover, the MLMC method runs multiple scenarios on different spatial and time meshes and then estimates the mean value of the mass fraction. The parallelization is performed in both the physical space and stochastic space. To solve every deterministic scenario, we run the parallel multigrid solver ug4 in a blackbox fashion. We use the solution obtained from the quasiMonte Carlo method as a reference solution.

MultipleRelaxation Runge Kutta Methods for Conservative Dynamical Systems(arXiv, 20230210) [Preprint]We generalize the idea of relaxation time stepping methods in order to preserve multiple nonlinear conserved quantities of a dynamical system by projecting along directions defined by multiple time stepping algorithms. Similar to the directional projection method of Calvo et. al., we use embedded RungeKutta methods to facilitate this in a computationally efficient manner. Proof of the accuracy of the modified RK methods and the existence of valid relaxation parameters are given, under some restrictions. Among other examples, we apply this technique to ImplicitExplicit RungeKutta time integration for the Kortewegde Vries equation and investigate the feasibility and effect of conserving multiple invariants for multisoliton solutions.