Nonlinear Solver Research Articles

High-order allpass cascade-connected phase-system design utilizing fractional polynomial errors

ABSTRACT This paper proposes a nonlinear optimization method for designing a high-order recursive allpass digital phase system that comprises several cascade-connected second-order (SO) allpass sections. This method first splits a given phase design specification (approximation target) into a set of equal sub-specifications, and then the equal sub-specification acts as the ideal phase for each of the SO sections. Since an SO section can be designed via employing a linear programming to minimize a linear fractional error function, this splitting process generates an excellent starting point (good initial values for the coefficients of the cascaded SO sections) for further minimizing the overall approximation error. This paper first derives a wide variety of fractional polynomial error functions for designing cascade-connected allpass phase systems of different high orders, and then shows that a high-order cascade-connected allpass phase system can be designed by utilizing a nonlinear solver to minimize the peak deviation of such a fractional polynomial error. Using the derived fractional polynomial error functions in the system design enables the system designer to find an excellent starting point for further minimizing the approximation error and thus reach a convergent design solution. This paper first derives a set of fractional polynomial errors, and then adopts the eighth-order cascade-connected phase system design to exemplify the effectiveness of the fractional-polynomial-error-based design tactic. Furthermore, the stability of the cascade-connected high-order allpass system can be readily checked via confirming that all of the cascaded SO sections have their poles inside the unit circle in the complex plane.

A novel LSTM-integrated non-intrusive ROM for reliability analysis of hysteretic systems with large stochastic dimension

Failure probability estimation of hysteretic systems subjected to random process excitation appears in many practical engineering applications, such as structures subjected to earthquake and vibration control. This task is computationally challenging due to the need for repeated solutions of a high-resolution model. This cost grows significantly with the stochastic dimension — characterized by the number of random variables, in this case, arising from the discretization of the random excitation. To address this issue, a novel non-intrusive reduced order model (ROM) is proposed in this paper. Accordingly, the complexity due to spatial discretization is first reduced using proper orthogonal decomposition. However, interpolation among response time histories in the resulting reduced solution space is computationally infeasible due to large stochastic dimensionality. An existing ROM developed by the authors is now improved by adopting a long short-term memory (LSTM) network for this interpolation. Furthermore, this adoption required two additional modifications in the mathematical structure of the ROM related to data compression and nonlinear coupling. The proposed methodology can be seamlessly integrated with black-box nonlinear solvers. Through detailed numerical studies on a beam on Winkler foundation and a multi-storied building frame subjected to stationary and non-stationary excitations, the proposed ROM is found to be accurate, efficient, and robust with respect to the frequency band of interest. The ROM has shown reduction of two orders of magnitude in the computational time. Furthermore, the proposed ROM is inherently parallelizable and holds promise to be scaled to other types of nonlinear systems with large stochastic dimensionality.

Simultaneous feeder reconfiguration and distributed generation planning in the presence of voltage-dependent and variable loads

In distribution systems, where distribution losses and the output of distributed generators (DGs) are significantly impacted by load power, an effective approach for reducing power loss is required, which the hybrid operation of DG units and grid reconfiguration can serve as a best alternative. Load power exhibits variability, altering alongside the voltage fluctuations which occur over time. Furthermore, the correlation between power demand and voltage relies on the type of load. Nonetheless, the incorporation of these important concerns into research on reconfiguration and distributed generation planning is rare. Only a limited number of papers have taken into account the voltage dependence and the type of time-varying loads in their respective models. Nevertheless, they proposed models with significant nonlinearity, requiring computation through nonlinear solvers or metaheuristic algorithms. Meanwhile, these models require the use of intensive linearization techniques to facilitate their implementation through linear solvers. High computational time is demanded by nonlinear solvers, while metaheuristic algorithms cannot guarantee the attainment of optimal solutions. Hence, the accurate modeling of load behavior holds importance in the active distribution system reconfiguration. In this paper, a proficient reconfiguration model is presented, which is both straightforward for implementation in conventional optimization tools and adept at identifying appropriate solutions for the reconfiguration and DG planning problem.

On the efficient non-linear solver for hydraulic fracturing and well cementing simulations based on Anderson acceleration

The aim of this study is to create a fast and stable iterative technique for numerical solution of a quasi-linear elliptic pressure equation. We developed a modified version of the Anderson acceleration (AA) algorithm to fixed-point (FP) iteration method. It computes the approximation to the solutions at each iteration based on the history of vectors in extended space, which includes the vector of unknowns, the discrete form of the operator, and the equation's right-hand side. Several constraints are applied to AA algorithm, including a limitation of the time step variation during the iteration process, which allows switching to the base FP iterations to maintain convergence. Compared to the base FP algorithm, the improved version of the AA algorithm enables a reliable and rapid convergence of the iterative solution for the quasi-linear elliptic pressure equation describing the flow of particle-laden yield-stress fluids in a narrow channel during hydraulic fracturing, a key technology for stimulating hydrocarbon-bearing reservoirs. In particular, the proposed AA algorithm allows for faster computations and resolution of unyielding zones in hydraulic fractures that cannot be calculated using the FP algorithm. The quasi-linear elliptic pressure equation under consideration describes various physical processes, such as the displacement of fluids with viscoplastic rheology in a narrow cylindrical annulus during well cementing, the displacement of cross-linked gel in a proppant pack filling hydraulic fractures during the early stage of well production (fracture flowback), and multiphase filtration in a rock formation. We estimate computational complexity of the developed algorithm as compared to Jacobian-based algorithms and show that the performance of the former one is higher in modelling of flows of viscoplastic fluids. We believe that the developed algorithm is a useful numerical tool that can be implemented in commercial simulators to obtain fast and converged solutions to the non-linear problems described above.

On the benefits of a multiscale domain decomposition method to model-order reduction for frictional contact problems

In this paper, the efficiency of a multiscale strategy based on a domain decomposition method (DDM) for model-order reduction of time-dependent frictional contact problems is presented. The proposed strategy relies on the LArge Time INcrement (LATIN) nonlinear solver combined with model reduction based on the Proper Generalized Decomposition (PGD). The LATIN presents a robust treatment of contact conditions, sharing similarities with augmented Lagrangian approaches, and naturally leads to a mixed DDM. In addition, the global space–time formulation of the method allows PGD-based model reduction to be used during computations, creating and enriching on-the-fly reduced bases per substructure to better track sliding fronts and propagative phenomena. The introduction of a multiscale strategy in the LATIN framework is consistent with the physics of contact problems, in which phenomena with different wavelengths interact: local solutions at contact interfaces presents high gradient effects with a short wavelength compared to the characteristic length of the structure. By taking advantage of this, the coarse scale problem of the strategy enables to capture efficiently the behavior of the problem at the structural level, focusing then on capturing the local contact variations at the contact interfaces. The most important features of the approach are shown comprehensively on a simple one-dimensional frictional contact problem. Then, its robustness and effectiveness are tested on a two-dimensional multibody frictional contact problem with more complex loadings. Guidelines are also given for choosing the parameters of the strategy, in particular those driving the construction of the reduced basis.

Solving nonlinear structural problems by associating line-search and path-following techniques

The demand for computational tools to simulate the nonlinear behavior of structures has intensified. Regarding nonlinear static or dynamic analysis, it is fundamental to use numerical strategies to trace the structure equilibrium paths, overcoming critical points (limit and bifurcations points). The nonlinear solvers must have a high level of efficiency in the two phases of the solution process (predictor and corrector), for each load step. In solving the nonlinear equations, it is quite common that Newton-Raphson's iterations do not converge or require an excessive number of iterations near equilibrium path critical points. Therefore, the line-search optimization technique appears as an additional sophistication. Basically, this technique aims to stagger the corrective displacements vector in the iterative phase, seeking to guarantee and accelerate the convergence of the process. This paper aims to verify the efficiency of the line-search technique coupled with Newton-Raphson iterations and different path-following methods and to verify their influence on the efficiency of the nonlinear solver. The effectiveness of the implemented line-search algorithm is verified by solving two slender structures with accentuated geometric nonlinearity. Such a resource is perceived to be triggered near load limit points (with more success when applied to structures with these critical points), accelerate the iterative process and increase the chances of convergence.

Jiezi: An open-source Python software for simulating quantum transport based on non-equilibrium Green's function formalism

We present a Python-based open-source library named Jiezi, which provides the means of simulating the electronic transport properties of nanoscaled devices on the atomistic level. The key feature of Jiezi lies in its core algorithm, i.e., self-consistent orchestration between the non-equilibrium Green's function (NEGF) method and a Poisson's equation solver. Beyond the construction of the tight-binding (TB) Hamiltonian with empirical parameters for conventional materials, the package offers a comprehensive framework for constructing the Wannier-based Hamiltonian matrix, enabling the investigation of novel materials and their heterostructures. To expedite the solution of NEGF systems, a methodology based on renormalization theory is proposed for reducing the dimension of the Hamiltonian matrix. Additionally, we adopt a non-linear Poisson equation solver with no analytical approximation in this software. The software facilitates seamless integration with external tools for geometry and mesh generation and post-processing. In this paper, we present the main capabilities and workflow by demonstrating with a simulation for the carbon nanotube field-effect transistor (CNTFET). Program summaryProgram Title: JieziCPC Library link to program files:https://doi.org/10.17632/nk79kbtww4.1Developer's repository link:https://github.com/Jiezi-negf/JieziLicensing provisions: GPLv3Programming language: PythonNature of problem: Simulates the quantum transport property of nano-scaled transistors based on the predefined device structure and the material composition.Solution method: Solves the coupled Schrödinger equation and Poisson equation by NEGF and finite element method.

Deep learning of first-order nonlinear hyperbolic conservation law solvers

In this contribution, we study the numerical approximation of scalar conservation laws by computational optimization of the numerical flux function in a first-order finite volume method. The cell-averaging inherent to finite volume schemes presents a challenge in the design of such numerical flux functions toward achieving high solution accuracy. Dense neural networks, as well as (piecewise) polynomials, serve as function classes for the numerical flux. Using a parametrization of such a function class, we optimize the numerical flux with respect to its solution accuracy on a set of Riemann problems for which we have closed-form solutions. By design, our resulting first-order finite volume schemes are conservative.We show empirically that the proposed schemes can achieve higher accuracy than the Godunov and Enquist-Osher method on a wide range of discretizations.The proposed method can be applied to the inverse problem, i.e., the data-driven synthesis of macroscopic physical laws for real-world flux phenomena.

Initial load factor adjustment through genetic algorithm for the generalized displacement control method: Implementation on non-rigid origami analysis

The nonlinear analysis of structural systems relies on effective nonlinear solution schemes, known as nonlinear solvers or path-following methods. This paper explores the challenge of determining the initial load factor in such analyses, with a particular focus on complex origami structures. This study presents a novel approach that leverages genetic algorithms (GA) to efficiently determine the initial load factor for the generalized displacement method. Building upon the concept of stiffness parameters, our method offers a systematic solution. In this work, GA serves as a tool to discover the appropriate initial load factor based on changes in the stiffness parameter, providing a more rational and user-friendly approach. To demonstrate the practical applicability of our approach, we apply it to a range of complex origami problems, including the folding of Miura-ori, folding of an egg-box origami, bending of an egg-box origami, and folding of an arc Miura-ori. These examples, along with a truss example, encompass highly nonlinear behaviors, including limit points, extreme stiffness fluctuations, snap-through, and snap-back. Our results showcase the effectiveness of our method in accurately tuning the initial load factor to capture and track intricate nonlinear behaviors. Furthermore, it offers the flexibility to fine-tune the distribution of points along the path and manage the resolution of the depicted equilibrium path. This versatility underscores its potential applications in fields such as deployable structures and complex large deformation structural analysis.

Mathematical Models for the Single-Channel and Multi-Channel PMU Allocation Problem and Their Solution Algorithms

Phasor measurement units (PMUs) are deployed at power grid nodes around the transmission grid, determining precise power system monitoring conditions. In real life, it is not realistic to place a PMU at every power grid node; thus, the lowest PMU number is optimally selected for the full observation of the entire network. In this study, the PMU placement model is reconsidered, taking into account single- and multi-capacity placement models rather than the well-studied PMU placement model with an unrestricted number of channels. A restricted number of channels per monitoring device is used, instead of supposing that a PMU is able to observe all incident buses through the transmission connectivity lines. The optimization models are declared closely to the power dominating set and minimum edge cover problem in graph theory. These discrete optimization problems are directly related with the minimum set covering problem. Initially, the allocation model is declared as a constrained mixed-integer linear program implemented by mathematical and stochastic algorithms. Then, the 0/1 integer linear problem is reformulated into a non-convex constraint program to find optimality. The mathematical models are solved either in binary form or in the continuous domain using specialized optimization libraries, and are all implemented in YALMIP software in conjunction with MATLAB. Mixed-integer linear solvers, nonlinear programming solvers, and heuristic algorithms are utilized in the aforementioned software packages to locate the global solution for each instance solved in this application, which considers the transformation of the existing power grids to smart grids.

Coupled 3D thermal-hydraulic code development for performance assessment of spent nuclear fuel disposal system

As a solution to the problem of spent nuclear fuels (SNFs), the disposal of SNF has gained attention from nations using nuclear energy because of hazards posed to the ecosystem. Among many proposed solutions, the most promising method is to dispose of SNF in a deep geological repository (DGR) which utilizes the multi-barrier concept developed by Finland and Sweden. Here, a new fully-coupled Thermal-Hydraulic (TH) code HADES (High-level rAdionuclide Disposal Evaluation Simulator) is developed using the MOOSE framework. This new code suggests basic numerical tools, such as a non-linear solver and finite element discretization, to assess the safety performance of disposal systems. The new TH code considered various TH behavior using Richards' flow approach, assuming gas pressure is constant. The HADES showed promising results when it was compared to various TH codes validated from DECOVAELX-THMC projects. When the single-canister model was utilized to estimate the TH behavior of the Korean Reference disposal System, although it showed significant saturation reduction due to the evaporation of water, the temperature was maintained under the thermal criteria limit, which is 100 °C. In addition, the new code estimated temperature and degree of saturation of the multi-canisters model, considering two or three canisters, it showed a slightly lower temperature, 5 °C, than the single-canister model. From these results, the following are concluded: (1) the new TH code contribute to an additional integrity by estimating TH behavior of KRS; (2) however, due to limitations in single-canister simulation, it is recommended to use multi-canisters simulation to estimate TH behavior accurately. Therefore, this model is anticipated not only to help licensing applications and estimation of various multi-physics phenomena and multi-canister at the disposal site.

An MFC-based friction damper with adjustable normal force: conception, modelling, and experiment

For classic friction dampers, the normal force on the contact interface is designed to be constant and cannot be easily adjusted post-manufacturing. The damping performance will attenuate when the working condition deviates from the target one. To overcome this issue, semi-active friction dampers are developed, where the normal force can be regulated by actuators to match the current working condition. However, the additional mass brought by the actuators limits their application prospects in the aerospace field. In this work, we propose a macro-fiber-composite-based adjustable friction damper (MFC-AFD) for typical thin-walled structures in aircrafts, which has the advantages of being lightweight, easy to install, and inherently fail-safe. The idea is to use the embedded MFC to drive the flexible structure of the damper to deform, thereby regulating the normal force. By adjusting the constant component of the normal force, the effective working range is expected to be enlarged, especially under non-target working condition. By controlling the second harmonic component, the optimal damping performance is expected to be further improved. We conduct the research using a combination of numerical simulations and experiments. To accurately and efficiently predict the damping performance of the MFC-AFD, a homogenization modeling method of MFC is proposed, and the corresponding experimental calibration technique is given. The harmonic balance-based nonlinear solver specifically for calculating the steady-response of electromechanically coupled structures is developed. Concerning the constant normal force control, a step-back to zero control law is established to eliminate the hysteresis behavior of MFC. In terms of the time-varying normal force control, a closed-loop control strategy in the frequency domain for the second harmonic component is proposed and realized. The vibration reduction effects of the MFC-AFD adopting the constant and time-varying normal force control laws are validated experimentally through a cantilever beam. Experiment results are in line with our expectations. When the vibration level is nearly 10 times of the target one, the reduction ratio recovers to 75% (i.e. close to the target working condition) by using the constant normal force control. There is an additional 32% vibration reduction compared to the optimal constant normal force when the time-varying normal force is properly controlled.

Boussinesq-type equations of hydroelastic waves in shallow water

Accurate computation of hydroelastic waves in shallow water is critical because many hydroelastic wave applications are nearshores, such as sea-ice and floating infrastructures. In this paper, Boussinesq assumptions for shallow water are employed to derive nonlinear Boussinesq-type equations of hydroelastic waves, in which non-uniform distribution of structural stiffness and varying water depth are considered rigorously. Application of Boussinesq assumptions enables complicated three-dimensional problems to be reduced and formulated on the two-dimensional horizontal plane, therefore the proposed Boussinesq-type models are straightforward and versatile for a wide range of hydroelastic wave applications. Two configurations, a floating plate and a submerged plate, are studied. The first-order linear governing equations are solved analytically with periodic conditions assuming constant depth and uniform stiffness, and the linear dispersion relations are subsequently derived for both configurations. For flexural-gravity waves of a floating plate, unique behaviours of flexural-gravity waves different from shallow-water waves are discussed, and a generalized solitary wave solution is investigated. A nonlinear numerical solver is developed, and nonlinear flexural-gravity waves are found to have smaller wavelength and celerity than their linear counterparts. For hydroelastic waves of a submerged plate, dual-mode analytical solutions are discovered for the first time. Numerical computation has demonstrated that a plate with decreasing submerged depth is able to transfer wave energy from the deeper water to the surface layer.

Efficient multimode vectorial nonlinear propagation solver beyond the weak guidance approximation

In this paper, we present an efficient numerical model able to solve the vectorial nonlinear pulse propagation equation in circularly symmetric multimode waveguides. The algorithm takes advantage of the conservation of total angular momentum of light upon propagation and takes into account the vectorial nature of the propagating modes, making it particularly relevant for studies in ring-core fibers. While conventional propagation solvers exhibit a computational complexity scaling as Nmode4, where Nmode is the number of considered modes, the present solver scales as Nmode3/2. As a first example, it is shown that orbital angular momentum modulation instability processes take place in ring-core fibers in realistic conditions. Finally, it is predicted that the modulation instability process is followed by the appearance of breather-like angular structures.

A highly parallel algorithm for simulating the elastodynamics of a patient-specific human heart with four chambers using a heterogeneous hyperelastic model

In this paper a highly parallel method is developed for simulating the elastodynamics of a four-chamber human heart with patient-specific geometry. The heterogeneous hyperelastic model is discretized by a finite element method in space and a fully implicit adaptive method in time, and the resulting nonlinear algebraic systems are solved by a scalable domain decomposition algorithm. The deformations of the cardiac muscles are quite complex due to the realistic geometry, the heterogeneous hyperelasticity of the cardiac tissue, and the myocardial fibers with active stresses. Moreover, the deformations in different chambers and at different phases of the cardiac cycle are very different. To simulate all the muscle movements including the atrial diastole, the atrial systole, the isovolumic contraction, the ventricular ejection, the isovolumic relaxation, and the ventricular filling, the temporal-spatial mesh needs to be sufficiently fine, but not too fine so that the overall computing time is manageable, we introduce a baseline mesh in space and a two-level time stepping strategy including a uniform baseline time step size to obtain the desired time accuracy and an adaptive time stepping method within a baseline time step to guarantee the convergence of the nonlinear solver. Through numerical experiments, we investigate the performance of the proposed method with respect to the material coefficients, the fiber orientations, as well as the mesh sizes and the time step sizes. For an unstructured tetrahedral mesh with more than 200 million degree of freedoms, the method scales well for up to 16,384 processor cores for all steps of an entire cardiac cycle.

Index‐aware learning of circuits

AbstractElectrical circuits are present in a variety of technologies, making their design an important part of computer aided engineering. The growing number of parameters that affect the final design leads to a need for new approaches to quantify their impact. Machine learning may play a key role in this regard; however, current approaches often make suboptimal use of existing knowledge about the system at hand. In terms of circuits, their description via modified nodal analysis is well‐understood. This particular formulation leads to systems of differential‐algebraic equations (DAEs), which bring with them a number of peculiarities, for example, hidden constraints that the solution needs to fulfill. We use the recently introduced dissection index that can decouple a given system of DAEs into ordinary differential equations, only depending on differential variables, and purely algebraic equations, that describe the relations between differential and algebraic variables. The idea is to then only learn the differential variables and reconstruct the algebraic ones using the relations from the decoupling. This approach guarantees that the algebraic constraints are fulfilled up to the accuracy of the nonlinear system solver, and it may also reduce the learning effort as only the differential variables need to be learned.

I-FENN with Temporal Convolutional Networks: Expediting the load-history analysis of non-local gradient damage propagation

In this paper, we demonstrate for the first time how the Integrated Finite Element Neural Network (I-FENN) framework, previously proposed by the authors [1,2], can efficiently simulate the entire loading history of non-local gradient damage propagation. To achieve this goal, we first adopt a Temporal Convolutional Network (TCN) as the neural network of choice to capture the history-dependent evolution of the non-local strain in a coarsely meshed domain. The quality of the network predictions governs the computational performance of I-FENN, and therefore we perform an extended investigation aimed at enhancing them. We explore a data-driven vs. physics-informed TCN setup to arrive at an optimum network training, evaluating the network based on a coherent set of relevant performance metrics. We address the crucial issue of training a physics-informed network with input data that span vastly different length scales by proposing a systematic way of input normalization and output un-normalization. We then integrate the trained TCN within the nonlinear iterative FEM solver and apply I-FENN to simulate the damage propagation analysis. I-FENN is always applied in mesh idealizations different from the one used for the TCN training, showcasing the framework’s ability to be used at progressively refined mesh resolutions. We illustrate several cases that I-FENN completes the simulation using either a modified or a full Newton–Raphson scheme, and we showcase its computational savings compared to both the classical monolithic and staggered FEM solvers. We underline that we satisfy very strict convergence criteria for every increment across the entire simulation, providing clear evidence of the robustness and accuracy of I-FENN. All the code and data used in this work will be made publicly available upon publication of the article.

Hybrid methods for solving structural geometric nonlinear dynamic equations: Implementation of fifth-order iterative procedures within composite time integration methods

This paper studies a class of hybrid methods that implement multi-point iterative procedures as the nonlinear solver within optimized composite implicit methods. These multi-point methods are employed to enhance convergence order and accuracy without requiring higher-order derivatives for solving structural geometric nonlinear dynamic equations. The promising approach provides an opportunity to delve deeper into its implications and explore the potential applications of multi-point techniques and composite methods across the analysis of dynamical systems. This study involves the utilization of two-point and five-point iterative methods with fifth-order convergence within three-step and four-step optimized composite implicit time integration methods. Moreover, a comprehensive hybrid scheme is introduced by integrating a family of optimized composite implicit methods with single-point and multi-point iterative methods. Then, the performance of the hybrid scheme is analyzed through structural examples, considering factors such as time step size, tolerance threshold, and the degrees of freedom of structures. The findings and numerical results illustrate that multi-point methods outperform the single-point method in terms of the number of iterations. Simultaneously, the four-step time integration method exhibits superior performance compared to the three-step time integration method, albeit with a marginal difference. It is observed that the hybrid method, consisting of the four-step time integration method and the five-point iterative method, performs better than other methods in terms of the number of iterations and errors for solving geometric nonlinear dynamic equations, especially in large structures. Thus, it stands out as a favorable choice for practical analyses and a possible candidate for generalization to other types of nonlinear problems in computational mechanics.

Solving a real world non-convex Quadratic Assignment Problem

The India-Japan Lighting company operates three plants in India. These three plants manufacture headlamps and taillights for the automotive industry. This study examines the facility location problem in one of these plants where 12 facilities must be placed in a two-column multi-row cellular layout. The machining sequences for the 20 parts conveyed among the 12 facilities were specified. The Quadratic Assignment Problem (QAP) is classified as an NP-hard problem in large instances. We modelled the specific instance as a QAP and reported the solution obtained by an easily available generalised reduced gradient (GRG) nonlinear solver and the solution obtained from the Gurobi optimiser. The Gurobi optimiser provides an excellent incumbent solution in quick time, but takes exponential time to reduce the duality gap.

Minimizing the Formation of Thermal Residual Stresses of the Elastic Zone in Functionally Graded Materials

Abstract: Functionally graded materials (FGM’s) are advanced composites consisting of two or more materials whose composition and structure changes gradually over the volume, leading to gradual change in its property. Typically one component is a ceramic, while the other is metal or metal alloy. In common embodiments, the composition can change gradually from all ceramic on one side to all metal on the other side, the ceramic contributing high resistance to temperature with the metal contributing high ductility. It’s a way to incorporate favorable properties of two materials into different locations of a single structure. In ideal cases the composition changes gradually but for manufacturing reasons, the changes occurs in homogeneous layers. The boundaries between the layers cause thermal residual stresses to arise, especially during material processing and in cases of cyclic loading. These thermal residual stresses are a result of differences in the thermal expansion coefficient of the varying compositions of materials. Thermal residual stresses may lead to cracking and ultimate failure of FGM’s. This study investigates the formation of thermal residual stresses for a Nickel-Alumina Ni-Al₂O₃ FGM. The Ni-Al₂O₃ system is chosen because it is one of the most common systems used in practice. This work explores the formation, impact and minimization of thermal residual stresses for a number of practical conditions that may arise during processing. All relevant properties are calculated using the rule of mixture equations. Abaqus, a non-linear FEA solver, is used for all the simulation work. In addition, the study highlights the need for examining elasto-plastic analysis to illuminate the process of crack development which is our future research goal.