Adjoint Formulation Research Articles

Inverse design of polaritonic devices

Polaritons, arising from the strong coupling between excitons and photons within microcavities, hold promise for optoelectronic and all-optical devices. They have found applications in various domains, including low-threshold lasers and quantum information processing. To realize complex functionalities, non-intuitive designs for polaritonic devices are required. In this contribution, we use finite-difference time-domain simulations of the dissipative Gross–Pitaevskii equation, written in a differentiable manner, and combine it with an adjoint formulation. Such a method allows us to use topology optimization to engineer the potential landscape experienced by polariton condensates to tailor its characteristics on demand. The potential directly translates to a blueprint for a functional device, and various fabrication and optical control techniques can experimentally realize it. We inverse-design a selection of polaritonic devices, i.e., a structure that spatially shapes the polaritons into a flat-top distribution, a metalens that focuses a polariton, and a nonlinearly activated isolator. The functionalities are preserved when employing realistic fabrication constraints such as minimum feature size and discretization of the potential. Our results demonstrate the utility of inverse design techniques for polaritonic devices, providing a stepping stone toward future research in optimizing systems with complex light–matter interactions.

Aerostructural Wing Optimization Using a Structural Surrogate in a Coupled Adjoint Formulation

As the number of disciplines included in the multidisciplinary design-optimization process continues to increase, it is envisioned that some of the disciplinary tools will take the form of surrogate models, whereas others remain physics-based, depending on the requirements and stage of the design process. To simulate this in the context of an aerostructural optimization of an aircraft wing, the work presented herein features a high-fidelity aerodynamic flow solver, while a surrogate is employed to model the wing structure. This approach includes the evaluation of the sensitivities of both the aerodynamic and structural disciplines, using a coupled-adjoint formulation to enable gradient-based optimization. An important aspect of the method is that the surrogate is trained only once, prior to the optimization, and held fixed throughout. The surrogate in effect parameterizes the structural design process, and outputs the weight and stiffness of an optimized structure, given inputs of geometry parameters and sizing loads. To minimize the number of surrogate inputs and enable the representation of the entire structural design space, parameterized loads are used to build the surrogate. The method is applied to the optimization of the NASA Common Research Model, illustrating the effectiveness of the new approach.

Inverse Scheme to Locally Determine Nonlinear Magnetic Material Properties: Numerical Case Study

We are interested in the determination of the local nonlinear magnetic material behaviour in electrical steel sheets due to cutting and punching effects. For this purpose, the inverse problem has to be solved, where the objective function, which penalises the difference between the measured and the simulated magnetic flux density, has to be minimised under a constraint defined according to the corresponding partial differential equation model. We use the adjoint method to efficiently obtain the gradients of the objective function with respect to the material parameters. The optimisation algorithm is low-memory Broyden–Fletcher–Goldfarb–Shanno (BFGS), the forward and adjoint formulations are solved using the finite element (FE) method and the ill-posedness is handled via Tikhonov regularisation, in combination with the discrepancy principle. Realistic numerical case studies show promising results.

Differentiable solver for time-dependent deformation problems with contact

We introduce a general differentiable solver for time-dependent deformation problems with contact and friction. Our approach uses a finite element discretization with a high-order time integrator coupled with the recently proposed incremental potential contact method for handling contact and friction forces to solve ODE- and PDE-constrained optimization problems on scenes with complex geometry. It supports static and dynamic problems and differentiation with respect to all physical parameters involved in the physical problem description, which include shape, material parameters, friction parameters, and initial conditions. Our analytically derived adjoint formulation is efficient, with a small overhead (typically less than 10% for nonlinear problems) over the forward simulation, and shares many similarities with the forward problem, allowing the reuse of large parts of existing forward simulator code. We implement our approach on top of the open-source PolyFEM library and demonstrate the applicability of our solver to shape design, initial condition optimization, and material estimation on both simulated results and physical validations.

Converted wave tomography based on inverse level set and adjoint formulation

SUMMARY Shear wave velocity (Vs) is a fundamental property of elastic media whose estimation from PS converted waves is challenging and requires modelling the boundary where P to S conversion occurs. This paper presents a PS tomography where seismic wave conversion/reflection points correspond to reflectors modelled with the level-set function set to zero [ϕ(x, z) = 0]. The proposed method aims for stable Vs inversion in a seismic acquisition setting using multicomponent receivers. Synthetic models simulating true Vs, Vp and the location of the geological reflector are used in the study. The inversion starts by locating a flat reflector, ϕ(x, z) = 0, which defines the zone Ω1 between the surface and the reflector, where the initial Vs and Vp fields are also set. To calculate the traveltimes of incident PT (P wave that propagates in Ω1 from source to the reflector), converted PS and reflected PP waves, for both observed and modelled data (forward problem), the methodology proposed by Rawlinson and Sambridge is adopted. This method uses the arrival times of the P waves, Tpt, from the seismic source at each reflector point as secondary sources generating the times Tps and Tpp. These times are calculated as a solution to the eikonal equation by using the Fast Marching method. The PS and PP residual times are minimized by updating Vs, Vp and ϕ(x, z) = 0 through adjoint variables designed from a formulation using Lagrange Multipliers in a variational context. The performance of the algorithm is evaluated for models with synclinal, sinusoidal and monoclinal reflector geometries using numerical tests considering the inversion of: (1) ϕ, given the true values of Vs and Vp; (2) ϕ and Vs, given the true value of Vp; (3) ϕ and Vp, given the true value of Vs and (4) the three parameters ϕ, Vs and Vp, simultaneously. Good results are obtained by inverting Vs and ϕ, given the true value of Vp. The simultaneous inversion of the three parameters exhibits promising results, despite the illumination problems caused by the different distribution of the PS, PP and PT time gradients due to the geometry of the reflectors and the acquisition setting (sources–receivers in the same plane). The proposed tomography estimates Vs and reflector positions which could help in statics corrections and improve the lithological characterization of near surface.

On the utilization of the adjoint method in microwave tomography.

In microwave imaging, the adjoint method is widely used for the efficient calculation of the update direction, which is then used to update the unknown model parameter. However, the utilization and the formulation of the adjoint method differ significantly depending on the imaging scenario and the applied optimization algorithm. Because of the problem-specific nature of the adjoint formulations, the dissimilarities between the adjoint calculations may be overlooked. Here, we have classified the adjoint method formulations into two groups: the direct and indirect methods. The direct method involves calculating the derivative of the cost function, whereas, in the indirect method, the derivative of the predicted data is calculated. In this review, the direct and indirect adjoint methods are presented, compared, and discussed. The formulations are explicitly derived using the two-dimensional wave equation in frequency and time domains. Finite-difference time-domain simulations are conducted to show the different uses of the adjoint methods for both single source-multiple receiver, and multiple transceiver scenarios. This study demonstrated that an appropriate adjoint method selection is significant to achieve improved computational efficiency for the applied optimization algorithm.

Optimal control of multilateral waterflooding wells in carbonate reservoirs with uncertainty consideration

The benefits of adjoint formulation in studying and improving the economic and production viability of reservoir waterflooding problems has been reported in several researches. However, there is little or no recorded report on its use for multilateral wells in uncertain waterflood reservoir conditions. This study aims to investigate the effect of geological uncertainty to multilateral well performance using adjoint formulation in oil reservoir waterflooding. First, an oil reservoir model with uncertainty assumption in permeability, Corey value, and physical dimension in reservoir compartmentalization was developed. Second, an optimal control sequence for the multilateral injection well was developed to optimize the net present value (NPV). The study revealed that with the obtained control sequence, reservoir compartmentalization notably faults in rock structure, significantly influenced the production capacity of multilateral wells which resulted in a higher NPV and oil rate when compared to other uncertainty cases.

A duality-preserving adjoint method for segregated Navier–Stokes solvers

Adjoint methods efficiently compute gradients for systems with many inputs and have been widely used for large-scale gradient-based optimization in fluid mechanics. To ensure optimization's numerical robustness, we need to develop an adjoint solution algorithm that has a similar, if not the same, convergence rate as the primal flow solver at each optimization iteration. This consistent primal-adjoint convergence behavior is also called duality-preserving (DP). Existing DP adjoint methods are mostly for fully coupled Navier–Stokes (NS) solvers with compressible flows. However, few DP adjoints have been proposed for segregated NS solvers, e.g., the semi-implicit method for pressure-linked equations (SIMPLE) algorithm solvers widely used for incompressible flow simulations. This study will fill this gap by deriving a fixed-point segregated adjoint formulation that fully preserves the convergence behavior of primal solvers. We first rewrite the steady-state segregated NS primal solution process into a fully coupled left-preconditioned Richardson format. Then, we transpose the preconditioner matrix to obtain a DP iterative adjoint formulation. In addition to algebraic derivation, we create a new graph representation that significantly streamlines the DP adjoint development for new segregated solvers. We prove that the proposed graph representation is equivalent to the above preconditioner transpose approach and guarantees duality. We evaluate our proposed algorithm using three cases with increasing complexity: a miniature-sized problem, an airfoil, and a wing. To quantify duality, we compare the eigenvalues between the primal and adjoint solvers and observe excellent agreements. We also evaluate to what extent various numerical settings (e.g., inner iteration tolerances and non-dual preconditioner) impact the eigenvalue distributions. Our proposed adjoint achieves machine-precision accurate gradient computation with competitive speed. Finally, we incorporate our adjoint solver into a large-scale gradient-based optimization framework and demonstrate its capability for wing aerodynamic shape optimization. The proposed fixed-point adjoint approach has the potential to make the adjoint solution more robust and efficient for any segregated NS solvers.

An Adjoint Method for High-Resolution EPMA Based on the Spherical Harmonics (PN) Model of Electron Transport

Electron Probe Microanalysis (EPMA) is a nondestructive technique to determine the chemical composition of material samples in the micro- to nanometer range. Based on intensity measurements of x-radiation, the reconstruction of the material composition and structure poses an inverse problem. The reconstruction methods currently applied are based on models that assume a homogeneous or a layered structure of the studied material. To increase the spatial resolution of reconstruction in EPMA the combination of a more sophisticated reconstruction method, which is based on a model that allows complex material structure, together with multiple measurements with varying beam configurations is required. We present a deterministic k-ratio model that is based on the PN model, an approximation of the radiative transfer equation for electron transport. Our goal is to approximate a maximum likelihood solution of the inverse problem using gradient-based optimization. We detail the application of the model in the context of algorithmic differentiation, in particular by deriving its continuous adjoint formulation. Algorithmic differentiation provides the flexibility to adapt the reconstruction method to various material parametrizations and thus to regularize and take into account prior knowledge. Through examples, we verify our implementation and demonstrate the flexibility of the reconstruction/differentiation method.

Adjoint-based control of three dimensional Stokes droplets

We develop a continuous adjoint formulation and implementation for controlling the deformation of clean, neutrally buoyant droplets in Stokes flow through farfield velocity boundary conditions. The focus is on dynamics where surface tension plays an important role through the Young-Laplace law. To perform the optimization, we require access to first-order gradient information, which we obtain from the linearized sensitivity equations and their corresponding adjoint by applying shape calculus to the space-time tube formed by the interface evolution. We show that the adjoint evolution equation can be efficiently expressed through a scalar adjoint transverse field. The optimal control problem is discretized by high-order boundary integral methods using Quadrature by Expansion coupled with a spherical harmonic representation of the droplet surface geometry. We show the accuracy and stability of the scheme on several tracking-type control problems.

Full-Waveform Inversion of Fiber-Optic VSP Data From Deviated Wells

Distributed acoustic sensing (DAS) technology permits measuring the strain or strain rate along a fiber deployed in a well due to seismic wave propagation. In vertical seismic profiling (VSP) acquisitions where seismic waves are excited at the surface and measured by subsurface receivers in the well, the DAS serves as a fast alternative to a conventional geophone array. Like conventional multicomponent VSP velocity measurements, DAS VSP data contain information about the subsurface properties that are encoded in seismic events and attributes such as compressional and shear direct arrivals, reflections, converted modes, multiples, and amplitude and phase variations. These seismic data can be used to obtain elastic properties of the formation and characterize the reservoir. Full-waveform inversion (FWI) extracts the elastic properties of the formation from the seismic wavefield by minimizing the misfit between real measurements and synthetically modeled wavefields. The synthetically generated wavefields are obtained by solving partial differential equations that accurately model the physical phenomena of the real waveforms. For conventional geophone sensors, FWI is formulated for the particle velocity field. In this paper, we present a formulation of the modeling and inversion specifically for DAS measurements as an averaged strain along the fiber. The algorithm does not require converting the data to velocities. The gradient of the misfit function, by adjoint formulation, is a cross correlation in time of the forward-propagated wavefield and backward-propagated residuals injected from the receivers. For conventional sensors, the residuals are injected as force terms, but for DAS data, the residuals are averaged in space first and then injected as moment tensor sources with radiation patterns determined by the well deviation. For multi-offset in-plane and out-of-plane VSP acquisitions, a two-dimensional (2D) algorithm has been developed that inverts for a 2D distribution of elastic medium properties and includes three-dimensional (3D) well deviation effects. This paper will present the results of applying the inversion to real multi-offset, highly out-of-plane VSP fiber data acquired in a deviated well. Simulations confirm that the well trajectory has a considerable impact on the data amplitude and must be included in the modeling and inversion to reproduce amplitude variations. The inverted compressional and shear velocities agree well with the reference model based on sonic data. Another real data example is a multi-offset walk-above acquisition where a targeted domain is located below a deviated portion of the well. The inversion matches reflections and produces an image in the target area.

Aerodynamic Shape Optimisation Using Parametric CAD and Discrete Adjoint

This paper presents an optimisation framework based on an open-source CAD system and CFD solver. In this work, the high-fidelity flow solutions and surface sensitivities are obtained using the primal and discrete adjoint formulations of SU2. This paper shows the direct use of CAD models for optimisation by developing a CAD system application programming interface and creating a link between the CAD-MESH-CFD analysis. A methodology to obtain geometric sensitivities is introduced, enabling the calculation of accurate gradients with respect to CAD variables and the deformation of the analysis mesh during the optimisation process. This methodology guarantees that the new surface mesh lies exactly on the CAD geometry. The optimisation framework is applied to a rectangular wing and a three section high-lift aerofoil configuration derived from the NASA CRM-HL configuration. Both geometries are created using FreeCAD. The performance objectives are to decrease the drag while constraining the lift to be above a desired value. The twist distribution of the wing was parameterised within the CAD system, allowing the minimisation of the induced drag by obtaining a nearly elliptical lift distribution. For the high-lift configuration, the position and rotation of the flap and slat were parameterised with respect to the original section; the final optimised positions yield a drag reduction of approximately 16.5%. These results show that the CAD parameterisation can be reliably used to obtain efficient optimums while operating directly on the CAD geometries.

Brachistochronous motion of a flat plate parallel to its surface immersed in a fluid

We determine the globally minimum time $T$ needed to translate a thin submerged flat plate a given distance parallel to its surface within a work budget. The Reynolds number for the flow is assumed to be large so that the drag on the plate arises from skin friction in a thin viscous boundary layer. The minimum is determined computationally using a steepest descent, where an adjoint formulation is used to compute the gradients. Because the equations governing fluid mechanics for this problem are nonlinear, multiple local minima could exist. Exploiting the quadratic nature of the objective function and the constraining differential equations, we derive and apply a ‘spectral condition’ to show the converged local optimum to be a global one. The condition states that the optimum is global if the Hessian of the Lagrangian in the state variables constructed using the converged adjoint field is positive semi-definite at every instance. The globally optimum kinematics of the plate starts from rest with speed $\propto t^{1/4}$ and comes to rest with speed $\propto (T-t)^{1/4}$ as a function of time $t$ . For distances much longer than the plate, the work-minimizing kinematics consists of an optimum startup, a constant-speed cruising, and an optimum stopping.

Efficient application of stochastic Discrete Well Affinity (DiWA) proxy model with adjoint gradients for production forecast

In this paper, we describe adjoint gradient formulation for the Operator-Based Linearization modeling approach. Adjoint gradients are implemented in Delft Advanced Research Terra Simulator (DARTS) framework and applied for history matching using a proxy methodology. Due to the application of adjoint gradients, the computational efficiency of the discrete well affinity (DiWA) proxy model for production forecast is significantly improved. That allows us to derive several important extensions. The proxy methodology is further extended and validated for 3D three-phase black-oil problems. The results show that the gradient-based regression can provide good history matching and reconstruct a true petrophysical characterization when the initial guess is generated based on highly reliable geological information. For cases with a limited or not sufficient geological characterization, an efficient stochastic application of DiWA proxy model is proposed. This approach consists of massive sampling procedures for collecting different realizations based on high-fidelity statistics with filtering. These realizations are generated stochastically because they are not conditioned to any production information but the basic geological statistics of the reservoir. The trained DiWA proxy model demonstrates a small deviation between the model response and the observation data. When applying the refined DiWA model for the training, the error between the model response and observation data can be further reduced. The forecast based on the trained model has slightly larger variability but the deviation is still reasonable. The enhanced DiWA methodology presents an efficient and robust technique for creating an ensemble of stochastic proxy models that can be used in production forecast, flow diagnostic, and optimization. • A new adjoint framework using an Operator-Based Linearization approach is proposed. • The framework is extended for the application of a three-phase 3D proxy model. • The characterization of the petrophysical information of reservoir is discussed. • An efficient sampling approach is applied to stochastic Discrete Well Affinity model.

Solute Transport Control at Channel Junctions Using Adjoint Sensitivity

Water quality control and the control of contaminant spill in water in particular are becoming a primary need today. Gradient descent sensitivity methods based on the adjoint formulation have proved to be encouraging techniques in this context for river and channel flows. Taking into account that most channels and rivers include junctions with other branches, the objective of this study is to explore the adjoint technique on a channel network to reconstruct the upstream boundary condition of the convection-reaction equation. For this purpose, the one-dimensional shallow water equations and the transport equation for a reactive solute are considered. The control is formulated through the gradient-descent technique supplied with a first-order iterative process. Both the physical and the adjoint equations are supplied with suitable internal boundary conditions at the junction and are numerically solved using a finite volume upwind scheme. The results reveal that the adjoint technique is capable of reconstructing the inlet solute concentration boundary condition in an acceptable number of iterations for both steady state and transient configurations using a downstream measurement location. It was also observed that the reconstruction of the boundary condition tends to be less effective the further away the measurement station is from the target.

Developing a Lung Model in the Age of COVID-19: A Digital Image Correlation and Inverse Finite Element Analysis Framework.

Pulmonary diseases, driven by pollution, industrial farming, vaping, and the infamous COVID-19 pandemic, lead morbidity and mortality rates worldwide. Computational biomechanical models can enhance predictive capabilities to understand fundamental lung physiology; however, such investigations are hindered by the lung’s complex and hierarchical structure, and the lack of mechanical experiments linking the load-bearing organ-level response to local behaviors. In this study we address these impedances by introducing a novel reduced-order surface model of the lung, combining the response of the intricate bronchial network, parenchymal tissue, and visceral pleura. The inverse finite element analysis (IFEA) framework is developed using 3-D digital image correlation (DIC) from experimentally measured non-contact strains and displacements from an ex-vivo porcine lung specimen for the first time. A custom-designed inflation device is employed to uniquely correlate the multiscale classical pressure-volume bulk breathing measures to local-level deformation topologies and principal expansion directions. Optimal material parameters are found by minimizing the error between experimental and simulation-based lung surface displacement values, using both classes of gradient-based and gradient-free optimization algorithms and by developing an adjoint formulation for efficiency. The heterogeneous and anisotropic characteristics of pulmonary breathing are represented using various hyperelastic continuum formulations to divulge compound material parameters and evaluate the best performing model. While accounting for tissue anisotropy with fibers assumed along medial-lateral direction did not benefit model calibration, allowing for regional material heterogeneity enabled accurate reconstruction of lung deformations when compared to the homogeneous model. The proof-of-concept framework established here can be readily applied to investigate the impact of assorted organ-level ventilation strategies on local pulmonary force and strain distributions, and to further explore how diseased states may alter the load-bearing material behavior of the lung. In the age of a respiratory pandemic, advancing our understanding of lung biomechanics is more pressing than ever before.

On tailoring fracture resistance of brittle structures: A level set interface-enriched topology optimization approach

We propose a fully immersed topology optimization procedure to design structures with tailored fracture resistance under linear elastic fracture mechanics assumptions for brittle materials. We use a level set function discretized by radial basis functions to represent the topology and the Interface-enriched Generalized Finite Element Method (IGFEM) to obtain an accurate structural response. The technique assumes that cracks can nucleate at right angles from the boundary, at the location of enriched nodes that are added to enhance the finite element approximation. Instead of performing multiple finite element analyses to evaluate the energy release rates (ERRs) of all potential cracks—a procedure that would be computationally intractable—we approximate them by means of topological derivatives after a single enriched finite element analysis of the uncracked domain. ERRs are then aggregated to construct the objective function, and the corresponding sensitivity formulation is derived analytically by means of an adjoint formulation. Several numerical examples demonstrate the technique’s ability to tailor fracture resistance, including the well-known benchmark L-shaped bracket and a multiple-loading optimization problem for obtaining a structure with fracture resistance anisotropy.

Discrete Adjoint-Based Aerodynamic Shape Optimization Framework for Natural Laminar Flows

An adjoint-based aerodynamic shape optimization framework for natural laminar flows is developed. The laminar to turbulent-transition onset is predicted by the correlation-based Bas–Cakmakcioglu transition model that is coupled with the Spalart–Allmaras turbulence model. A discrete adjoint implementation is subsequently developed. Automatic differentiation is utilized to construct the partial derivatives in the discrete adjoint formulation. The turbulence and transition models are fully coupled into the sensitivity derivative evaluations, as well as into the objective function evaluations. The sensitivity derivatives evaluated by the discrete adjoint solver are validated against those of the finite central difference method. The discrete adjoint-based aerodynamic shape optimization framework developed for natural laminar flows is successfully employed to optimize aerodynamic characteristics of the NLF(1)-0416 airfoil and a low-aspect-ratio wing.

The continuous adjoint cut-cell method for shape optimization in cavitating flows

A continuous adjoint formulation for shape optimization in steady-state, cavitating flows is developed. A Transport Equation-based mixture model, extended with the Kunz cavitation model, is implemented to incorporate phase transition due to cavitation and for which the adjoint equations are derived. Flow and adjoint equations are discretized on 2D Cartesian meshes with cut-cells. In the cut-cell method, Cartesian cells intersected by the solid boundaries are reshaped by discarding their solid part, resulting in cells with an arbitrary number of sides. Emphasis is laid on the accurate computation of geometric sensitivities at the cut-cells, required by the adjoint method to compute the objective gradient. The sensitivity derivatives obtained via the adjoint method are compared with finite differences for the purpose of validation. The proposed adjoint formulation is assessed by performing three shape optimizations with two different objectives, namely cavitation suppression and hydrofoil lift maximization over isolated hydrofoils.

Adjoint EM Sensitivity Analysis for Fast Frequency Sweep Using Matrix Padé via Lanczos Technique Based on Finite-Element Method

Sensitivity analysis is important for electromagnetic (EM)-based design. The existing adjoint EM sensitivity analysis methods have to solve large systems of EM equations repetitively for different frequencies. This article addresses this situation and proposes to speed up the EM sensitivity analysis over a frequency range by solving EM equations at only a single frequency. A new adjoint EM sensitivity analysis algorithm for the fast frequency sweep using the matrix Padé via Lanczos (MPVL) technique based on the finite-element method (FEM) is proposed in this article. MPVL is incorporated to relate the information of one frequency to the information of multiple frequencies. A large system of EM equations is then solved at a single frequency to predict the sensitivity information for the entire frequency band. Adjoint formulations are further derived to avoid the effect of the number of design variables. The adjoint EM sensitivity analysis using the proposed technique can obtain the same accuracy as the existing techniques while taking less time by avoiding repetitively solving large systems of EM equations for different frequencies and different design variables. The proposed technique is demonstrated by three EM examples of microwave components.