Structural optimisation of diffusion driven degradation processes

In this article, we propose an optimisation framework that can contribute to the prevention of premature failure or damage to building structures and can thereby strengthen their longevity. We concentrate on structures that are contaminated by chemical substances and that have strong destructive effects on the material. The aim of this contribution is a mathematical algorithm that allows the optimisation of a structure exposed to chemical influences and thus the assurance of the static load capacity. To achieve this, we present a coupled mechanical-diffusion-degradation approach embedded in a finite element (FE) framework. Furthermore, we integrate an optimisation algorithm to reduce material degradation. In this paper, we establish shape optimisation of a structure with a gradient based optimisation algorithm.


Introduction
Engineering structures are dimensioned according to standards. The maximum existing stresses are evaluated and the material load-bearing capacity is examined, whereby environmental influences are only marginally considered. However, negative influences can additionally change the material composition and thus the mechanical load-bearing capacity over time. Most changes in the internal structure of materials are associated with diffusion processes. One example is the long-term effect of calcium leaching in concrete, where "pure water creates concentration gradients which lead to the diffusion of Ca ions from the pore water and the subsequent degradation of underground concrete" (Choi and Yang 2013).
In the following, an approach is presented which allows the calculation of mechanical and chemical influences on structures. Furthermore, it is possible to calculate an Division of Solid Mechanics, Lund University, Lund, Sweden optimal geometry that reduces the damage of long-term effects caused by environmental influences. Therefore, this approach integrates the structural analysis into a structural optimisation algorithm. Within the structural analysis, we couple the diffusion of chemical substances with the mechanical behaviour of the material. We assume permeable structures which allow gradient-based diffusion of concentrations. In addition, we postulate that chemical substances have an influence on the material of the structure, since they can trigger material degradation. To ensure this, we focus on mechanical degradation processes corresponding to negative growth processes using multiplicative decomposition of the deformation gradient. Since this is a space-and time-dependent problem, we embed the material description in a finite element (FE) framework, by using an isoparametric concept for space discretisation and a Newmarkbeta approach for temporal discretisation. Furthermore, we develop an algorithm for structural optimisation. In detail, we present a shape optimisation with a gradient-based calculation, which contains information about the sensitivities of the parameters that influence the mechanical behaviour.
Using representative examples, we introduce the basic idea of the model and examine practical aspects where components are at risk of being exposed to chemical concentrations and yet their strength must be guaranteed.
In order to determine the approach for optimisation of diffusion driven degradation problems, the following topics are applied and embedded in a short overview of the state of the art. In the context of structural mechanics, the topics coupled problems and growth, respectively degradation processes, are investigated in more detail. For the embedding in optimisation, shape optimisation and sensitivity analysis are applied. Kuhl (2005) presents a general overview of coupled problems and the numerical implementation. Furthermore, he outlines in detail the coupling of diffusion processes to mechanical behaviour. A coupled theory of fluid permeation and large deformations for elastomeric materials, concentrating on a thermodynamically consistent derivation of the constitutive relations and the resulting partial differential equations, is presented by Chester and Anand (2010). Menzel and Kuhl (2012) summarise growth and remodelling models for living structures. Mechanical growth and remodelling can be modelled either with a constitutive approach, a kinematic approach or a combination of both. Growth processes can be described by evaluating the timedependent change in mass, density or volume of a structure. On the one hand, the constitutive approach concentrates on a thermodynamically consistent evaluation of the mass source and the mass flow, which enables the calculation of change in mass or density (Cowin and Hegedus 1976;Harrigan and Hamilton 1992;Epstein and Maugin 2000;Kuhl et al. 2003;Menzel 2005). The kinematic approach, on the other hand, allows the calculation of the variable mass or volume by applying a multiplicative decomposition of the deformation gradient (Rodriguez et al. 1994;Chen and Hoger 2000;Ambrosi and Mollica 2002;Gleason and Humphrey 2005;Menzel 2007). A useful alternative is the combination of the constitutive and kinematic approach, one example for this combination is presented by (Ganghoffer 2010), where the calculation of surface growth in biological tissue is presented. He uses the multiplicative decomposition of the deformation gradient, as originally presented by Rodriguez et al. (1994), and applies a thermodynamically consistent approach to establish an evolutionary law for growth velocity. The multiplicative decomposition of the deformation gradient leads to a "growth tensor describing the local addition of material and an elastic tensor characterizing the reorganization of the body", see Ganghoffer and Plotnikov (2014). Therefore, the development of the growth tensor, the so-called transplant tensor, is introduced as a state variable in the framework of finite elasticity.
In Gérard et al. (2002), a simplified numerical model of calcium leaching is presented, which concentrates on the chemical state variable calcium and the kinetics of the leaching process. This model considers time-dependent chemical processes and mentions the possibility of a simple coupling in a finite element (FE) algorithm. Kuhl (2005) presents a chemo-mechanical model for the numerical simulation of calcium leaching in concrete. He emphasises the effects of the chemical degradation on the pore volume and the mechanical stability of concrete. A damage function is used for constitutive modelling of the chemo-mechanically degraded material.
An overview on structural optimisation methods applied to discretised linear-elastic structures is given in Christensen and Klarbring (2008). Structural optimisation, which is applied to growth or degradation processes, is usually solved by evolutionary optimisation algorithms (ESO), see Simon (2013) for details. The evolution models are based on the concept of gradually removing inefficient material from a structure. The disadvantage of ESO concepts is that optimisation starts from one defined reference configuration and that the result quickly ends in local minima. Thus, it is not possible to find an optimal solution with a single ESO algorithm, see Vrugt and Robinson (2007). However, the application of shape optimisation to growth or degradation processes could help solve the problem, but has not yet been sufficiently investigated. Similar to the embedding of a growth approach in an optimisation algorithm, shape and topology optimisations are applied to damage models. For example, an optimal damage distribution is computed with a shape optimisation algorithm that embeds an isotropic gradient-enhanced damage model, see Guhr et al. (2020). Suresh et al. (2018) present a topology optimisation method which is connected to a fatigue model. The model enables the computation of optimised topologies, taking fatigue at high cycles as a limiting condition into account. Moreover, Noël et al. (2017) develop a level set-based topology optimisation framework to reduce damage in the context of structural design. Barthold (2002) presents the basic principles of design changes necessary for structural optimisation and their effects on the structural response on a continuum. He introduces the local convective approach containing local coordinates and derives it from a differentiable manifold. Furthermore, he emphasises the importance of this approach for obtaining information about the kinematic relation required for numerical methods such as the finite element method (FEM) and computer-aided geometric design (CAGD). A detailed discussion of efficient strategies for calculating sensitivity of design parameters is described in Barthold and Stein (1996) where, in addition, the method of variational sensitivity analysis is presented.
The above briefly presented literature motivates the structural analysis which contains a coupled mechanicaldiffusion-degradation approach. The special feature of this model is the combination of the coupled problem with a shape optimisation approach so that the failure limit of structures affected by mechanical loads and chemical impacts can be improved. By embedding the degradation model in the optimisation algorithm, this approach represents an effective alternative to evolutionary algorithms.

Continuum model
We present the continuum model including coupling between diffusion and mechanical behaviour. First we introduce the extended kinematic framework. We then formulate the partial differential equations which are determined by balance equations and constitutive equations together with boundary and initial conditions. The constitutive equations, including mechanical stresses, particle diffusion and chemical reactions, are determined by a thermodynamically consistent framework based on the evaluation of the entropy inequality.

Kinematics
A graphical illustration of the applied kinematic framework is presented in Fig. 1, which provides a multiplicative decomposition of the deformation gradient into an elastic and a growth part. We assume that the temporal evolution of the mass density takes place in the reference configuration. An infinitesimal volume element is represented by dV in the reference configuration, by dv d in the degradation space and by dv in the actual configuration. We introduce positions in the reference configuration X, in the actual configuration x and in the degradation space x d . Furthermore, ϕ t describes the mapping of the reference particles X onto their spatial position x = ϕ t (X, t) with t representing time. We introduce the convective tangent vectors G i , g i and h i . The contravariant basis vectors that determine the dual basis are defined by G i · G j = δ j i , g i · g j = δ j i and h i · h j = δ j i . With this, we introduce the deformation tensors which map infinitesimal line elements represented in different configurations, i.e. F, F d and F e deformation gradient : (1) Fig. 1 Graphical illustration of the kinematic concept These gradients are important two-point tensors that allow transformations between objects in relation to the respective configurations. We apply a multiplicative decomposition of the deformation gradient F into an elastic part F e and a degradation part F d as introduced in Himpel et al. (2005), Kuhl et al. (2004), and Menzel and Kuhl (2012). According to Lubarda and Hoger (2002), we use an isotropic approach for the degradation fraction of the deformation gradient with the stretch ratio ν, to be specific whereby 1 denotes the second order identity tensor and the mass densities ρ 0 and ρ * 0 respectively, being introduced as this work proceeds. Using the elastic part of the deformation gradient, we can introduce the elastic right Cauchy Green tensor as Moreover, the material time derivative of the elastic right Cauchy Green tensor is obtained bẏ witḣ F e = l F e and d = 1 2 (l T + l) , wherein l = gradẋ is the spatial velocity gradient. Figure 1 introduces the definition of the mass densities referred to the respective configurations. In the reference configuration, ρ * 0 represents the initial and ρ 0 the referential mass density, which arises as a result of material degradation. The mass density ρ d is referred to the degradation space and ρ t to the actual configuration. Mass exchange is realised by a mass sink term R 0 per unit volume in the reference configuration which triggers the degradation from the initial mass contribution dM and results in the mass contribution dm, i.e.
This leads to time-dependent update of the initial mass density ρ * 0 to the resulting referential mass density ρ 0 , to be specific Within the kinematic concept, we assume that the initial referential mass density and the mass density referred to the degradation space coincide, i.e.
The mapping between the configurations, see Fig. 1, leads to a change in volume, whereas mass is not influenced by such mappings, i.e. using (8) results in Moreover, J e = det F e > 0 denotes the determinant of the elastic part of the deformation gradient and J = det F > 0 as well as (2), follow by analogy.

Degradation approach
The degradation model is implemented by a combination of a kinematic and a constitutive approach. Thus, in Section 2.1, the multiplicative decomposition of the deformation gradient is performed and the degradation space is introduced. Within the constitutive approach, a mass exchange is established which leads to the change of referential mass density ρ 0 . In return, the evolving mass density is used to calculate the degradation part of the deformation gradient, see (2). We assume that the mass degradation is caused due to chemical concentrations c γ , which are defined by the corresponding molar density ρ γ and the molar mass M γ , i.e.
Thus, we apply the constitutive approach of the mass sink term as follows Assuming only one concentration to be related to the chemical interaction, we can summarise the degradation impact of the concentrations as This paper does not address chemical processes in detail.
For further details on a model with more complex chemical processes, the reader is referred to, e.g. Gérard et al. (2002). In the following, the concentrations are introduced as an additional degree of freedom in the continuum.

Balance equations
In order to implement the above-introduced problem numerically, the following section presents the necessary balance approaches to solve the coupled multi-field problem. On the one hand, we need the balance of mass to calculate the concentrations, and on the other hand, the momentum balance to calculate the displacement. The concentration distribution depends on spatial conditions and the stresses on mass density which, in turn, is influenced by concentrationrelated mass reduction. For this reason, this is a strongly coupled problem. Furthermore, the balance of energy and entropy is needed to determine the constitutive material model.

Balance of mass
We evaluate the balance of mass for the macroscopic body and for the chemical concentrations. The reduction of the mass in the reference configuration leads to a change of mass with dṁ = R 0 dV and{ •} denoting the material time derivative. Using the transformations between the configurations, we can represent the balance of mass on the actual configuration as wherein 'div' denotes the spatial divergence operator. The balance of mass of the chemical concentrations contains the material time derivative of the concentrationsċ γ and the flux of the concentrations j γ , i.e.

Balance of linear momentum
The balance of linear momentum is shown below for the actual configuration and contains the Cauchy stress tensor T. Although we apply mass exchange, we can neglect effects of the impulse resulting from the degrading mass, since the process of degradation is assumed to be very slow. This results in wherein additional volume forces as well as acceleration are neglected. From the balance equation of moment of momentum, we can derive that the Cauchy stress tensor is symmetric, so that T T = T.

Balance of energy
The balance of energy consists of the temporal derivative of internal energyĖ, and kinetic energyK, which correspond to mechanical energy change dW, thermal energy change dQ and chemical flow on the surface E c γ , resulting in the following equationṡ whereinė is the material time derivative of the volume specific internal energy, r is the volume specific heat source, q is the heat flux density, μ γ is the chemical potential and da = n da with spatial outward unit surface vector n. The time derivation of the kinetic energy has no impact on the balance of energy as we neglect acceleration and consider only slow progression of degradation. The surface part can be reformulated as follows, Using the mass balance of chemical concentrations from (14), the local form of the balance of energy reads

Balance of entropy
Entropy is summarised in the second law of thermodynamics and states that entropy never decreases in a closed system. Accordingly, entropy always increases or remains constant. In this context, the entropy inequality is adapted as whereṡ denotes the material time derivative of the volumespecific entropy and where is the absolute temperature. The surface part can be reformulated as follows, with div 1 q = q · grad 1 + 1 div q = − 1 2 q · grad + 1 div q .
Inserting the reformulation into the entropy inequality (20), we obtain the local form Under isothermal conditions, i.e. q = 0 as well as˙ = 0, and with the local form of the energy balance, see (19), we end up with the restriction

Constitutive formulation
We perform a Legendre(-Fenchel) transformation between the extensive thermodynamic entropy s and the conjugate quantity, the temperature . This results in the so-called Helmholtz energy ψ, i.e.
By introducing the Helmholtz free energy into the entropy inequality (24) and by selecting the elastic right Cauchy Green tensor and the concentrations to derive the energy process with ψ(C e , c γ ), the following results are obtained We can separate the approach for the Helmholtz function into a mechanical ψ M and a chemical ψ C part In order to ensure the entropy inequality, the reversible contributions are evaluated and the Cauchy stress T and the chemical potential μ γ are determined, i.e.
wherein (4) is used. From the irreversible part of the entropy inequality we motivate the flux of the concentrations where D ≥ 0 is introduced as the diffusion coefficient.

Specifications of energy contributions
The mechanical energy is described by adopting a hyperelastic Neo-Hooke material, ψ Neo , which contains material parameters μ and λ, to be specific with invariants I C e = tr C e and J e = √ det C e = √ III C e , which depend on the elastic right Cauchy Green deformation tensor C e . In addition, we specify the chemical part ψ C as where μ 0 γ is the constant standard potential and where R is the gas constant. Inserting the mechanical and chemical energy contributions into (28) results in with the spatial Karni-Reiner strain tensor K e , i.e.
Furthermore, the chemical potential follows as

Illustration of stress state
In the evaluations of the examples in Sections 4 and 6 we use the principal Cauchy stresses T 1 and T 2 in the two-dimensional space considered, namely with the invariants of the Cauchy stress I T and II T , i.e.

Finite element formulation
The previously presented continuum framework for the coupled mechanical-diffusion-degradation model leads to a set of coupled differential equations, which are timedependent and highly nonlinear. For further procedure, a numerical approximation, the FEM, is used. We apply the Galerkin method, whereby the balance equations are represented in their weak form and weighted by independent test functions. For the fully discrete problem, we solve for the displacements u = x − X and the concentrations c γ on the basis of the balance of linear momentum and the balance of mass for the concentrations.

Weak form
In the following we present the weak form of the balance of momentum weighted by the independent test functions for the displacements δu and the weak form of the balance of mass of the concentrations weighted by the independent test functions for the concentrations δc γ . The equations are posed in the reference configuration. On the one hand, the weak form of the balance of linear momentum results in wherein the mapping of the Cauchy stress T to the first Piola-Kirchhoff stress tensor P is applied with the transformation P = J T F -T and wherein N is the material normal unit vector on the Neumann boundary. On the other hand, the weak form of the balance of mass is given as with J F -T N dA = n da. We describe the Dirichlet boundary conditions as follows and the Neumann boundary conditions as with J γ = J j γ F -T . The initial conditions are given with u(t 0 ) = u 0 and c γ (t 0 ) = c 0 γ .

Discretisation in space
For space discretisation, we apply the isoparametric concept which is based on approximating geometry, displacement and concentrations by the same set of ansatz functions h I (ξ ). The discrete form of the test functions for the displacement δu h and for the concentrations δc h γ results in wherein NN denotes the number of nodes per element and where ξ represents the local coordinates. In this work, we apply an eight-noded element description, i.e. two-dimensional Serendipity elements under plane strain conditions. Moreover, approximations of the degrees of freedom as well as of all related gradient operations follow straightforwardly.

Discretisation in time
The simulation of the diffusion of the concentrations requires a discretisation in time for which we apply the Newmark-beta method. Within the considered time interval, we approximate a constant average acceleration of the concentrationsc β , i.e.
with the previous acceleration of the concentrationsc γ (t) and the acceleration to be approximated in the present time stepc γ (t + t). Based on this, one obtainṡ

Numerical analysis examples
In this section, we discuss the properties of the proposed diffusion controlled degradation model in the context of representative examples focusing on the basic coupling mechanism between diffusion impact and structural response.

Material behaviour
In the following, the correlations between mechanical and degradation processes are investigated on the basis of an analysis including one finite element and two different sets of boundary conditions as shown in Fig. 2. Furthermore, we apply the material and loading parameters shown in Table 1. The degradation process is identically established for both examples, at this point independently of chemical concentrations, with the change of the referential mass density ρ 0 , as shown in Fig. 3. With this at hand, the degradation part of the deformation gradient F d can be evaluated with (2).
For mechanical impact, the first example is uniaxially loaded with a displacement u = 0.5 cm in y-direction, wherein ε is the strain in y-direction. For comparison, the second example is loaded with a constant traction with resulting force F = 0.015 MN in y-direction, see Fig. 2b. The determinant of the degradation contribution to the deformation gradient, J d , follows directly from the timedependent approach for the change of referential density, see Fig. 4, so that we observe the same referential mass density evaluations for both boundary conditions. However, the total deformation changes differently in the two examples. For the boundary condition displayed in Fig. 2a, the total deformation represented by J is constant, because the structure does not allow shrinkage. The elastic deformation  The grey colour shows the initial state, the black colour the deformation state and the blue colour the degraded structure: displacement loading (left, a) and force loading (right, b) by analogy with Fig. 2 changes inversely so that, i.e. J e = [J d ] −1 J. The stress T yy increases, because the total deformation is fixed, whereas the degradation part decreases. In comparison to the first example with fixed Dirichlet boundary conditions, the example illustrated in Fig. 2b results in a degradation process and the total deformation represented by J decreases over time. In return, the elastic deformation J e and the stress T yy are constant. The related states of deformation are additionally illustrated in Fig. 5.

Diffusion driven degradation
In this example, we discuss the coupled effects of diffusion leading to a degradation of the material. This example is motivated by the idea to analyse chemical influence on hollow concrete blocks, since concrete is a porous medium that allows the inflow of concentrations and is susceptible to chemical degradation. The boundary value problem is illustrated in Fig. 7a, which shows a structure with a hole. The hole of the structure is defined by the parameters (s 1 , s 2 ), which represent the axes of an ellipse. We consider a concentration inflow from the left side of the structure. The choice of parameters, as stated in Table 2, is based on values common for concrete and Table 2 Material parameters for the structure with the hole = 50 cm s 1 = 10 cm s 2 = 5 cm the orders of magnitude for chemical diffusion processes, although no specific chemical process is described here. The material degradation triggered by molecular processes with concentrations occurs very slowly. Therefore, the calculation is accelerated by considering the diffusion rate per day (d) and a time-dependent increase in concentrations per day, as shown in Fig. 6. In total, the computation considers a time period of 4 d with a time step size of 1 d. A stable mechanical environment is enabled by fixed displacements as shown in Fig. 7a. The concentration inflow on the left side leads to the contour plot given in Fig. 7b. Figure 8a illustrates the degradation induced by the concentrations and in Fig. 8b the impact of the deformation on the first principal stress T 1 is evaluated. Since material degradation is triggered by a high concentration of  The maximum of the first principal stress is located in this area close to the hole. Overall, this example clearly shows the coupling between diffusion and deformation.

Structural optimisation framework
This section outlines the connection between structural analysis and optimisation framework. Furthermore, the information required for the optimisation problem, such as the objective function, the constraints and design parameters, are presented. The main programme runs in the numerical computing environment of MATLAB which contains a link to the open software gmsh to create a mesh on the one hand and an interface to a Fortran code using MEX-file interfaces on the other. The linearisation of the weak form, in discrete form, is implemented in a Fortran based FE code and the MEX-file transfers the data to the workspace. The assembly and calculation are implemented in MATLAB.

Sketch of algorithmic framework
An overview of the algorithm is shown in Fig. 9. The algorithm can be divided into three main sections, the boundary value problem, structural analysis and structural optimisation. First, the boundary value problem is relevant for the definition of the model problem by means of parametric geometry description. The information about the boundary value problem provides the input for the structural analysis. In addition, the solution of the structural analysis computes a structure deformed by chemical and mechanical loads, which is used as initial design in the structural optimisation. In this paper, continuum mechanical Fig. 9 Illustration of the algorithmic framework quantities such as the evaluation of the stress restriction and geometry parameters are selected for the objective function and the constraints. These variables are calculated within the framework of structural analysis and directly passed on to the mathematical solver. The structural optimisation minimises the defined objective function while maintaining the given constraints and provides new design parameters. The optimisation task is solved with the help of a MATLAB toolbox.
Furthermore, the technical implementation can be briefly summarised. The global programme runs in MATLAB, while the information at element level is implemented in a Fortran code and embedded via MEX interfaces. Information at element level retrieves the coupled equations for the description of the problem presented, the variation formulations with the discrete weak forms and the gradient information for the FEM solution. In addition, an interface to the open software gmsh, a finite element mesh generator, is generated during the framework of the boundary value problem. In this paper, geometrical values are chosen as design parameters; therefore, a parametric mesh design is implemented. The functions introduced as part of the optimisation process are described in the typographic style 'italics', e.g. objective function J , constraints g, design parameters s.

Optimisation problem
The structural optimisation problem is solved by using the nonlinear optimisation function fmincon provided by MATLAB toolbox, see (MathWorks 2019). This solver finds the minimum of a restricted nonlinear and multivariable function. The general optimisation problem follows with an objective function J (v, s), nonlinear inequality constraints g(v, s), upper and lower limit values s u and s l for a set of design parameters, i.e. min J (v, s) : g(v, s) wherein v ∈ {u, c γ } are the field variables and s geometric design parameters. In this paper, we deal with the two shape optimisation problems as follows: 1. Minimisation of the maximum first principal stress T max 1 by changing geometrical parameters, with the constraint of a maximum loss of area. 2. Minimisation of the area by changing geometrical parameters, with the constraint of an upper limit for the maximum first principal stress T max 1 . In both optimisation setups, structural analysis is applied to calculate the objective and constraint functions, since there is a dependency on the field variables. Further specifications on the objective function and constraint are provided in Sections 6.1 and 6.2. Fmincon uses the 'sqp-legacy' algorithm to solve the optimisation task, wherein the gradients of the objective function and constraints are obtained numerically by finite differences. The Hessian matrix is iteratively integrated using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) procedure. Further details can be found in (MathWorks 2019). The application of the finite difference method is a precise but time-consuming method. In the outlook, we address a more efficient approach for future work.

Numerical optimisation examples
In the following, two examples are presented for shape optimisation of structures that are loaded by chemical concentrations. The first example is taken from Section 4.2 and shows a hollow concrete block which is loaded by chemical substances. We calculate the optimised shape of the brick in such a way that the maximum stress caused by the concentrations is reduced while still retaining material. The second example is inspired by a mechanically loaded bridge, additionally loaded with chemical substances, which can be caused by, for example calcium leaching. The goal of this example is a structure with minimised material and limited stresses.

Optimisation of a structure with a hole
Using the example of the structural analysis in Section 4.2 with the material parameters from Table 2, an optimised shape is discussed below. We use the sum of the Gaussian point values of the first principal stress in the maximum loaded element T max 1 as the objective function J (v, s) and the change of the area as inequality constraints g(v, s) with wherein NG is the number of Gaussian points, A ini is the initial area, A is the actual area andĀ is a threshold. The optimisation problem follows as The axes of the ellipse (s 1 , s 2 ) are the design parameters. The optimisation algorithm applies the input parameters listed in Table 3.
The following diagramme, see Fig. 10, shows the iteration course of the optimisation to minimise the objective function while fulfilling the constraints. The solver requires a total of 9 iterations until an optimal shape is found to minimise the average first principal stress in the Gaussian points to 0.083 MN cm −2 .
After 18,468 s of computing time, the optimisation leads to the following new design parameters, i.e. (Fig. 12). Fig. 10 Iteration of the optimisation solver 'sqp-legacy', which shows the decrease of the objective function, i.e. the first principal stress where the length of the axes of the hole change to reduce the first principal stress. The contour plot of the first principal stress in Fig. 11 shows that the change of the hole parameters, in other words the design parameters, leads to an overall decrease of the first principal stress. Thereby, the constraint is fulfilled, i.e.
This example illustrates how shape has a major influence on the effects of degrading concentrations and that damage can be minimised by changing the shape. This observation is supported by Fig. 13 which shows the temporal course of the deformation based on the determinant of the degradation gradient, J d , at the maximum loaded point. The curve is compared on the one hand with the initial shape, and on the other hand with the optimised shape, showing that the optimised form leads to a lower degradation.

Remarks on numerical investigation
In this paper, the mathematical optimisation is solved utilising the MATLAB function fmincon, which calculates the minimum of a constrained nonlinear multivariable function, see MathWorks (2019). The user can chooses between different algorithms to solve the task. In this section, the influence of the solution algorithm on the optimised result is evaluated. For this purpose, the result of the example 'Optimisation of a structure with a hole' is evaluated using three different algorithms. The 'sqp-legacy' and 'active-set' algorithm are based on sequential quadratic programming (SQP) method. Thereby, SQP is derived using the Newton method and taking into account inequality constraints. The two algorithms differ in their implementation, e.g. they apply different definitions for the strict feasibility with respect to bounds or the choice of the solution algorithm for the subproblems. In contrast, the third algorithm 'interior-point' combines two different approaches to solve the optimisation task. The algorithm uses Newton steps or conjugated gradient steps depending on the solution of each iteration step. Table 4 compares the efficiency of the algorithms and shows that algorithm 'sqp-legacy' runs most efficiently, Furthermore, Fig. 14 shows the optimised design taking into account different mesh sizes. The identical optimised parameters with s= [5; 5.525] are determined, but the computing time increases with the refinement of the mesh. In detail, the finer mesh shown in Fig. 14 needs 3.3 times more calculation time than the coarse mesh.

Optimisation of a bridge-like structure
This example is inspired by a bridge with mechanical loads under environmental influences, such as calcium leaching. This scenario can occur, for example during longterm exposure to pure water, which triggers the diffusion of calcium ions. To illustrate whether environmental influences are taken into account or not, the calculation is listed in the first step without and in the second step with the influence of chemical concentrations. In order to reduce the material costs, the objective J (v, s) of the problem is to minimise the area within a plane strain setting while holding a threshold for the first principal stress, i.e. g(v, s). The optimisation problem is obtained as follows wherein the area is calculated by the sum of the total element volumes. We apply a threshold for the maximum first principal stress T max 1 , thus the total optimisation problem follows with, i.e. min J (v, s) : g(v, s) ≤ 0 constraints s l ≤ s ≤ s u limit values.
To save computing time, the symmetry of the structure is utilised and the calculation is performed on half of the system using symmetry boundary conditions on the right side, see Fig. 15. Table 5 presents the applied material parameters. The applied forces with F = 10 kN lead to the first principal stress distribution, which is displayed in the contour plot in Fig. 18b. The formation of a tensile area in red and compression area in blue becomes visible.
The lower edge of the structure is defined by a Bspline function with four control points, where the design parameters are the control points that allow vertical displacement, see Fig. 16, i.e.
Utilisation of symmetry must also be taken into account for the calculation of design parameters during optimisation. For this purpose, the condition x 3 = x 4 is introduced.
The parameters for the optimisation algorithm are shown in Table 6, where s u and s l are the matrices containing the minimum and maximum allowable change of design parameters. Figure 17 shows the iteration of the optimisation solver with the decrease of the objective function, i.e. the area of the structure. After 6 iterations the objective function converges, the side condition is fulfilled and a local minimum is given with the following design parameters s = 1; 0.893; −0.1372; −0.1372 .
The design parameters reveal that the largest saving occurs in the area of the least stress, i.e. in the area where neither compressive nor tensile stress is present, whereby more material is required in the middle of the structure to ensure load-bearing capacity. In total, the optimised design   Fig. 18 represents the first principal stress in the optimised and initial design. In a second step, the optimised bridge is additionally loaded by chemical concentrations. In this example, we focus on the general effect of any chemical concentrations that trigger material degradation. Additionally, to the mechanical force, two concentration inflows are located on the top of the structure and are provided by a concentration increase of 0.33 mol m −3 per time step, see Fig. 19.
As already mentioned in the example in Section 4.2, the chemically induced degradation process is accelerated by running the simulation in time steps of days. The contour plots in Fig. 20 show the distribution of the concentrations and the resulting influence on the determinant of growth after 9 days. Figure 21 shows the time progression of the first principal stress in one nodal point, P 1 , where the maximum first principal stress occurs. The first principal stress increases due to the increased concentrations. The previously defined maximum stress of T max 1 = 120 kN m −2 can no longer be maintained.
The aim of the optimisation is to save as much material as possible while still ensuring the load capacity. For this reason, smaller deviations in stress due to diffusion processes are also relevant, since otherwise the load-bearing capacity can no longer be guaranteed over a long load duration.  By analogy with the previous example, the points of the polygon chain at the lower edge of the structure are design parameters. Using the same parameters for the optimisation algorithm as provided in Table 6 and starting from the initial design as illustrated in Fig. 15, the optimisation yields the following optimal design parameters s = 0.933; 0.793; −0.196; −0.196 .
The iteration process illustrated in Fig. 22 shows how the design parameters reduce the objective function by 17.4 % while maintaining the constraints. The increased maximum first principal stress due to the degradation process changes the optimal design and allows for less material saving, as shown in Fig. 23. It is clear that the consideration of environmental conditions, such as calcium leaching, are necessary to predict long-term performance. Fig. 22 Iteration of the optimisation solver 'sqp-legacy' for the example coupled to chemical impact, which shows the decrease of the objective function, i.e. the area Fig. 23 Comparison of the results: the green coloured area shows the optimal result when chemical influence is considered and the red area shows the result without degradation processes

Summary
In this paper, a coupled mechanical-diffusion-degradation model is presented. The degradation process is derived by a multiplicative split of the deformation gradient and a constitutive approach for the development of growth, which assumes chemical concentrations as a trigger for degradation. The numerical FE framework is briefly presented. Furthermore, the embedding of a structural optimisation framework is outlined. The applicability of the model is presented for a structure with a hole and a beam for which a practical reference is outlined. The main focus of the examples lies in the analysis of the influence of longterm acting chemical concentrations, which can influence the mechanical stress and can keep certain upper limits of the stress by constructive changes.
For future work, an alternative to the numerical finite difference method will be applied, namely the variational sensitivity analysis as outlined in Barthold and Stein (1996), which accelerates the simulation time.
Overall, it is shown that the model can provide an optimal design taking into account long-term effects from concentrations that damage the material, while still maintaining certain limits for the load-bearing capacity. The algorithm offers the possibility to integrate different chemical processes, to calculate the interaction with the mechanical behaviour and to solve optimisation problems.
Funding Open Access funding enabled and organized by Projekt DEAL.

Replication of results
The described method is outlined in the essential steps and the continuous formulation is elaborated. Furthermore, the most important algorithmic features are discussed. The examples contain all necessary information about the material model and geometry parameters required for replication.

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