Computational modelling and data-driven homogenisation of knitted membranes

by   Sumudu Herath, et al.
University of Moratuwa

Knitting is an effective technique for producing complex three-dimensional surfaces owing to the inherent flexibility of interlooped yarns and recent advances in manufacturing providing better control of local stitch patterns. Fully yarn-level modelling of large-scale knitted membranes is not feasible. Therefore, we consider a two-scale homogenisation approach and model the membrane as a Kirchhoff-Love shell on the macroscale and as Euler-Bernoulli rods on the microscale. The governing equations for both the shell and the rod are discretised with cubic B-spline basis functions. The solution of the nonlinear microscale problem requires a significant amount of time due to the large deformations and the enforcement of contact constraints, rendering conventional online computational homogenisation approaches infeasible. To sidestep this problem, we use a pre-trained statistical Gaussian Process Regression (GPR) model to map the macroscale deformations to macroscale stresses. During the offline learning phase, the GPR model is trained by solving the microscale problem for a sufficiently rich set of deformation states obtained by either uniform or Sobol sampling. The trained GPR model encodes the nonlinearities and anisotropies present in the microscale and serves as a material model for the macroscale Kirchhoff-Love shell. After verifying and validating the different components of the proposed approach, we introduce several examples involving membranes subjected to tension and shear to demonstrate its versatility and good performance.



There are no comments yet.


page 15

page 16


Gaussian process modelling using UQLab

We introduce the Gaussian process modelling module of the UQLab software...

Deformation and failure in nanomaterials via a data driven modelling approach

A data driven computational model that accounts for more than two materi...

Gaussian Process (GP)-based Learning Control of Selective Laser Melting Process

Selective laser melting (SLM) is one of emerging processes for effective...

Tensor Basis Gaussian Process Models of Hyperelastic Materials

In this work, we develop Gaussian process regression (GPR) models of hyp...

Online tuning and light source control using a physics-informed Gaussian process Adi

Operating large-scale scientific facilities often requires fast tuning a...

Approximate Solution Approach and Performability Evaluation of Large Scale Beowulf Clusters

Beowulf clusters are very popular and deployed worldwide in support of s...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Knitting is one of the most efficient and widely used techniques for producing fabric membranes. The recent advances in computational knitting make it possible to produce large complex three-dimensional surfaces in one piece without seams [26, 37]. On a modern programmable flat-bed knitting machine it is possible to produce even non-developable surfaces by a local variation of the stitch pattern consisting of an interlooped yarn. Most promisingly, the yarn can be replaced or integrated with electroactive or conductive yarns to produce novel interactive textiles with sensing and/or actuation capabilities [24, 25, 34]. Knitted membranes are usually very flexible because the primary deformation mechanism for the yarn is bending rather than axial stretching. Their unique stretchability and drapability properties make knitted membranes appealing as a reinforcement in composite components [21, 14]. If needed, the stiffness can be increased by inserting straight high-strength fibres during the knitting process. Such reinforced membranes have been recently used as a formwork in architectural engineering [36]. Considering the recent advances in knitting, there is a need for efficient computational approaches for the analysis of large-scale knitted membranes.

For large-scale analysis of knitted membranes, the computational homogenisation approaches which take into account the deformation of the interlooped yarn on the microscale are crucial. There is an extensive amount of literature on finite element-based computational homogenisation of heterogeneous solids, see e.g. the reviews [33, 13, 39]. Two-scale homogenisation often referred to as FE, combined with a yarn-level and a membrane-level finite element model has been also applied to woven and knitted fabrics [29, 12, 22, 11]

. The boundary conditions of the microscale representative volume element (RVE) are given by the membrane deformation and in turn the averaged yarn stress in the RVE yields the membrane stress. However, such schemes are inefficient for knitted membranes with large deformations because of the need to solve a nonlinear problem at each quadrature point of the membrane. It is increasingly apparent that for nonlinear problems two-scale homogenisation must be considered in combination with a data-driven machine learning model 


. The model can be trained in an off-line learning phase by solving the microscale problem for a sufficiently rich set of deformation states. Subsequently, the trained model provides a closed-form constitutive equation that is used in the macroscale model. As a machine learning model, for instance, neural networks or GPR models have been used 

[2, 20, 31, 51]. Alternatively, it is possible to formulate the macroscale finite element problem directly on the training data set bypassing the need to train a machine learning model [17, 35]

. GPR regression is a well-studied Bayesian statistical method and provides as such a principled approach to dealing with epistemic and aleatoric uncertainties and issues such as overfitting 

[38]. Therefore, we approximate the response of the microscale yarn model with GPR. Their limitation to relatively small dimensions and number of training points is not relevant for this paper.

Owing to the relative slenderness of the yarn on the microscale and the membrane on the macroscale, they are best modelled as a rod and a shell, respectively. Evidently, this leads in comparison to a 3D solid model to an immense reduction in the number of degrees of freedom. The geometrically exact rod theories pioneered by Simo et al. 

[43, 42] provide a consistent and efficient framework for modelling rods undergoing finite deformations. Later contributions to geometrically exact rod theories include [3, 19, 16]. While in most of these classical rod models the transverse shear deformations are taken into account, in more recent Euler-Bernoulli type models they are neglected [27, 1]. The omission of the transverse shear effectively sidesteps the shear locking problem which is a major impediment in the analysis of slender rods. These new models are usually discretised using smooth B-spline basis functions due to the presence of the higher-order displacement derivatives in the energy functional. The Euler-Bernoulli type rod model introduced in this paper takes into account the stretching, bending and torsion of the yarn as well as the non-frictional rod-to-rod contact between yarns. Similar to the yarn, the fabric membrane can be modelled with geometrically exact shell theories going back again to Simo et al. [40, 41]. As for rods, in more recent Kirchhoff-Love type shell models the transverse shear deformations are neglected and the weak form is discretised using smooth B-splines or subdivision basis functions [7, 8, 18]. Although the bending resistance of knitted membranes is usually very small, a suitably chosen shell bending energy term can serve as a regularisation for wrinkling under compressive stresses [8, 5]. In this paper we make use of the Kirchhoff-Love subdivision shell implementation introduced in [6, 23].

The outline of the paper is as follows. In Section 2 we introduce our finite deformation rod model, its discretisation with B-splines as well as the treatment of rod-to-rod contact. Subsequently, in Section 3 we discuss the proposed microscale yarn-level RVE model for computational homogenisation. This model is verified and validated with experimental and numerical results from literature. In Section 4, we introduce the data-driven GPR model and describe its training with the yarn-level RVE model. Finally, in Section 5, we first discuss the training of the GPR model and then analyse membranes subjected to tension to assess the accuracy of the obtained data-driven constitutive model.

2 Finite deformation analysis of yarns

In this section we summarise the governing equations for the finite deformation Euler-Bernoulli rod model for the yarn. The presented equations are without loss of generality restricted to rods with circular cross-sections. We take into account rod-to-rod contact by enforcing the non-penetration constant with the Lagrange multiplier method.

2.1 Kinematics

The geometry of the rod is described by a set of circular cross-sections connected by their line of centroids. In accordance with the Euler-Bernoulli assumption, equivalent to the Kirchhoff-Love assumption for shells and plates, transverse shear is neglected so that cross-sections remain always normal to the line of centroids.

Figure 1: Geometric description of a spatial rod in its reference (left) and deformed (right) configurations. The two configurations are defined using the line of centroids  and 

and the respective covariant basis vectors.

The position vectors of material points in the reference and deformed configurations and are parametrised in terms of the convective coordinates  as


where  and  denote the lines of centroids parameterised by , see Figure 1. In turn, the cross-section with the radius  is parameterised by  and . The tangent vectors to the line of centroids are given by


The two orthonormal directors  and  are chosen so that they satisfy in the reference configuration


Here and in the following Greek indices take the values  and summation over repeated indices is assumed. The two orthonormal directors  in the deformed configuration are obtained by rotating the reference configuration directors with a rotation matrix ,


This, in combination with the Euler-Bernoulli assumption, ensures that the two directors  and  in the deformed configuration satisfy


The rotation is composed of two rotations,


That is, the reference directors are mapped to the deformed directors in two steps using an intermediary configuration with . The matrix describes a rotation by an angle  about the unit tangent and is according to the Rodrigues formula given by


with the identity matrix

and the skew-symmetric matrix


Subsequently, we use the smallest rotation formula [27, 9] for the second matrix , which maps  to the unit tangent vector ,


To derive the strains corresponding to the assumed kinematics (1

), we consider the Green-Lagrange strain tensor of 3D elasticity


Here and in the following Latin indices take the values . The covariant basis vectors  and  and the contravariant basis vectors  and are defined as


where  is the Kronecker delta. The corresponding two metric tensors  and  are given by


After introducing the assumed kinematics (1) in (10) and some algebraic simplifications we obtain for the components of the strain tensor


We identify as the membrane, bending about , bending about and torsional shear strains, respectively. Moreover, due to the Euler-Bernoulli assumption the strain components and are zero and the in-plane shear strains and are induced only by torsional shear strain .

2.2 Equilibrium equations in weak form

The potential energy of a rod with the line of centroids  and the cross-section  occupying the volume  in its reference configuration takes the form


where is the strain energy density and  is the potential of the externally applied forces. At equilibrium the potential energy of the rod is stationary, i.e.,

with the external virtual work  and the internal virtual work

where is the second Piola-Kirchhoff stress tensor. The strain tensor of the rod (13) depends on the displacement of the line of centroids


and the rotation angle . Hence, we can write for the internal virtual work (15b) more succinctly


where and are the internal forces conjugate to the virtual displacements  and rotations .

As a material model we use the isotropic St Venant-Kirchhoff model with the strain energy density


and the fourth-order constitutive tensor


where and are the two Lamé parameters [4]. The contravariant metric tensor  is determined from the relation .

To derive analytical expressions for the internal forces we first introduce the rod strain tensor (13) and the constitutive equation (19) in the internal virtual work (17). Subsequently, we integrate over the rod cross-section analytically to obtain the internal forces


Here, the axial force

, the bending moments

and the torque are defined as


where , and are the cross-section area, second moment of area and torsional constant of the circular rod, and is its Young’s modulus. The tedious but straightforward derivation of internal forces and their derivatives are summarised in Appendix A.

2.3 Finite element discretisation

We follow the isogeometric analysis paradigm and use univariate cubic B-splines to discretise the lines of centroids and in the reference and deformed configurations. We choose smooth B-splines because the bending strains in (13) require at least continuous smooth basis functions. The kinematic relationship (16) is restated after discretisation as


where the B-spline basis and their coefficients correspond to the control vertices on the discretised rod centreline.

The discretised weak form of the equilibrium equations (15) yields after linearisation an algebraic system of equations, which we solve with the Newton-Raphson scheme. For linearisation the gradients of the internal forces (20) are required, see Appendix A.2.

2.4 Yarn-to-yarn contact

Figure 2: Contact between two rods with the lines of centroids  and . Between the points  and  the distance is minimum and has the value .

We use the Lagrange multiplier method to consider pointwise non-frictional contact between two circular rods. By closely following Wriggers et al. [50] and Weeger et al. [48], we add to the total potential energy in (14) the contact potential energy


where the Lagrange multiplier represents the repulsive normal force between the two rods, and the non-positive gap function  depends on the minimum distance between the two rods. The distance between the lines of centroids of the two rods and is given by


and its minimum by


This minimum can be determined using the Newton-Raphson scheme. Ultimately, the non-positive gap function  is defined as


where is the radius of the rods.

The first variation of the contact potential energy (23) gives its contribution to the weak form of the equilibrium equations (15). This contribution takes the form


For further details on our contact implementation we refer to Herath [15].

2.5 Verification of the rod model

Figure 5: One-sided clamped helicoidal spring with tip loading. Problem description (a) and deformed shapes under different loads (b).
Figure 6: Comparison of the obtained tip displacements of the helicoidal spring with Bauer et al. [1].

We consider a membrane-bending-torsion interaction problem to verify the accuracy of the presented rod formulation. A spatial helicoidal spring is clamped at one end and a vertical load of magnitude 25 kN is applied at the other end, see Figure (a)a. The reference geometry of the line of the centroids is given by


The material and geometric properties of the spring including the reference director orientations  are given in Figure (a)a. The tip load is increased from zero to 25 kN in 10 uniform load steps. In Figure (b)b the deflected shapes for five different load levels are shown. As can be seen in Figure 6 the obtained tip displacements are in excellent agreement with the results presented in Bauer et al. [1], in which a slightly different approach was used to parameterise the rotations of the two rod directors .

3 Computational Homogenisation

3.1 Microscale analysis

A characteristic representative volume element (RVE) as depicted in Figure 7 is chosen to represent the periodic microstructure of a weft-knitted membrane. The intricate spatial arrangement of the yarns in the RVE is modelled, as in recent works [10, 49], using the approximate geometry proposed by Vassiliadis et al. [45], see Appendix B. For alternative geometric descriptions see [47] and references therein. First the yarn centrelines are defined and then the yarns are generated by sweeping a circle along those lines.

Figure 7: Characteristic weft-knitted RVE selection and the definition of fibre geometric and material parameters.

On the macroscale we model the membrane with subdivision shell finite elements [8, 6, 23] and consider for homogenisation only the in-plane membrane response. The very small out-of-plane bending stiffness contribution of the rods to the membrane bending stiffness is not taken into account.

In first-order homogenisation, the macroscale membrane deformation gradient  is used to define the boundary conditions on the RVE. Here and in the following the macroscale and microscale quantities are denoted by subscripts and , respectively. It is assumed that the volume averaged deformation gradient of the rod in the RVE and the membrane deformation gradient  are equal, that is,


where is the RVE volume and  is the rod volume within the RVE in the reference configuration. Without going into details, the deformation gradient of the line of centroids is given by the kinematic assumption (1). The Hill-Mandel lemma states that the macroscale and microscale work for an RVE must be equal


where  is the macroscale first Piola-Kirchhoff membrane stress and  is the microscale first Piola-Kirchhoff rod stress. The integral represents the internal virtual work of the rod within the RVE and is given by (17). For the lemma (30) to hold only certain types of RVE boundary conditions can be chosen, including Dirichlet and periodic; see [13, 52] for details. Moreover, as a consequence of this lemma the macroscopic first Piola-Kirchhoff stress  can be obtained from the volume averaged internal energy density


where the set  denotes the deformation gradients satisfying the RVE boundary conditions. The microscopic energy density  of a rod with an isotropic St Venant-Kirchhoff material is given by (18).

Figure 10: Prescribed periodic boundary conditions and the respective deformed configuration of an RVE subjected to biaxial tension (a) and shear (b). At the six boundary finite element nodes, the four labels describe the boundary conditions for the three displacements and the one rotation, i.e. , where refers to the applied displacements in the direction. The labels and denote periodic and zero displacement or twist constraints, respectively.

In this work we choose as RVE boundary conditions the periodic boundary conditions depicted in Figure 10. Due to the orthotropy of the RVE response, the biaxial stretching and shear response are decoupled. Therefore, we consider the two cases separately and choose for each different boundary conditions. Boundary displacements in the thickness direction are zero to simulate plane stress conditions.

3.2 Verification and validation of the microscale model

To verify and validate our microscale model we consider RVEs under wale-wise and coarse-wise uniaxial tension and shear as reported in Dinh et al. [10] and Weeger et al. [49]. For a detailed problem description we refer to [10] and [49]. In our computations the yarn is discretised with 128 rod finite elements. In Dinh et al. [10] the yarn is modelled as a 3D solid and in Weeger et al. [49] as a discretised Euler-Bernoulli rod with the collocation method. As evident from Figure 11 our results are in excellent agreement with the experimental and computational results obtained with other approaches. Throughout this paper, we use the yarn geometry and material parameters given in the above-mentioned two papers.

Figure 11: Comparison of the obtained homogenised stresses for an RVE under course-wise and wale-wise uniaxial tension and pure shear with the experimental and numerical results reported in Weeger et al. [49] and Dinh et al. [10]. The unit conversion for stresses is .

4 Data-driven homogenisation

4.1 Review of Gaussian process regression

Gaussian process regression is a statistical inference method rooted in Bayesian statistical learning. In the following, we briefly review the key steps in Gaussian process regression. For further details we refer to Rasmussen et al. [38]. To begin with, we assume as a prior in the Bayesian formulation a random function given by the zero-mean Gaussian process


The respective covariance function is chosen to be the (stationary and isotropic) squared exponential function


where and

are the yet to be determined scaling and characteristic lengthscale hyperparameters. For convenience, we define

as the vector of hyperparameters. The infinitely smooth covariance function (33) encodes our prior assumptions about  before observing any training data. Other choices are possible, see [38, Chapter 4].

Next, we consider a database comprised of a known training dataset and testing points. Each data point is given by an input vector and a corresponding scalar observation , with . All training data points are collected in (, ) and all testing data points with unknown in (, ). Moreover, we define the covariance matrix with the components . The discretisation of the Gaussian process (32

) for the considered training and test data is given by the multivariate Gaussian distribution


The probability density of the unknown data

at the given locations conditioned on the known training data reads


Hence, the best estimate for

is given by the expectation


and the uncertainty of the estimate is given by the covariance matrix


The covariance matrix  of size  is dense for the covariance kernel (33). If the number of training points  become too large to invert , alternative approaches leading to sparse covariance matrices need to be considered [38, Chapter 8].

To estimate the hyperparameters , we consider the marginal likelihood, or the evidence,


with the likelihood and the prior , where are function values evaluated at , c.f. (32). The marginal likelihood is the probability of observing the data for given hyperparameters . Hence, it suggests itself to choose the hyperparameters so that the probability observing the given training data  is maximised. In practice, it is numerically more stable to compute the maximum of the log marginal likelihood given by


Note that the is a monotonic function so that


We use the scikit-learn Python library [32] to find  with a gradient-based optimisation algorithm. The marginal likelihood is a non-convex function and certain care has to be taken to find the global maximum. After the optimal hyperparameters are determined they are used in computing the mean (36) and covariance (37).

4.2 Gaussian process homogenisation

Figure 12: Data-driven computational homogenisation and material design framework using Gaussian process regression.

The implemented data-driven homogenisation framework is illustrated in Figure 12 and closely follows Bessa et al. [2]. Data-driven homogenisation begins with constructing a response database for the RVE. This response database is designed by defining the microscale boundary value problem (BVP) using three types of design variables, namely geometric properties , material properties and boundary conditions . In this paper, we concentrate on the design variable relevant to homogenisation. It is straightforward to include and which lead to the material design of knitted membranes.

Next, we determine the hyperparameters of the GPR model using -fold cross-validation rather than directly using the  obtained with (40) considering the entire training dataset. The GPR model is iteratively evaluated by using folds for training and the remaining one for testing. We choose which provides a good compromise between accuracy and efficiency. As error metrics for the predicted strain energy and stress resultants we use the correlation of determination and the mean squared error (MSE) given by


where is the number of the test points,

is the target output,

is the predicted output and is the expectation of the target output given by the GPR model. When evaluating the MSE of the predicted stress tensor, we take the square root of the squared sum of the MSE of each stress resultant


The GPR model training is implemented in a Python environment using scikit-learn machine learning library [32]. For each test fold, we determine the corresponding hyperparameters by maximising the log marginal likelihood (40). It is expected that the obtained hyperparameters have similar values for all test folds. In practice, this is not the case because the log marginal likelihood is usually a non-convex function. Therefore, we vary during optimisation the upper and lower limits of the hyperparameters so that all test folds yield similar hyperparameter values. Once this is achieved, we record the hyperparameters  of the GPR model, which has the lowest MSE, and save the model for subsequent macroscale analyses.

Lastly, we integrate our in-house thin-shell solver [7, 23] with the trained GPR model to simulate the homogenised knitted membrane. As shown in Figure 12, for given macroscale deformations (strains) the trained GPR model yields the corresponding predicted stress resultants and their tangents to the macroscale shell solver.

4.3 Macroscale analysis

For plane-stress membrane deformation, a response database is created with the three in-plane components of the Green-Lagrange strain  as the inputs and RVE volume averaged strain energy as the output, i.e. and  in Section 4.2. With a slight abuse of notation  is here, in contrast to Section 2.1, a two-dimensional second order tensor. In the four-dimensional response database, points have the coordinates . The strain components and correspond to a biaxial deformation state that is decoupled from the shear deformation state with , see Section 3.1. The RVE strain energy  is obtained by considering a biaxial strain state with  and and a shear deformation state with  and summing up their strain energies. We use the validated strain limits of Figure 11 in constructing the response database. Thus, allowing for a 5 compressive strains, we define the design variables of the response database as,


We solve the RVE problem with the boundary conditions stipulated by the strain states defined in (43) and store the strain components and respective volume averaged energies in the response database. After passing the GPR model training and testing phase, as shown in Figure 12, the strain energy density  of the plane-stress deformation is predicted for a given deformed state of the RVE. Plane stresses and stress tangents are computed by differentiating (36), that is,


5 Examples

5.1 GPR model training and testing

We start with presenting the error metrics of GPR model training and testing for varied sizes of response databases. GPR model training errors are observed to take values and . Hence, we present only the error metrics of testing datasets. All the following computations are performed on the Intel Core i5-4590 CPU @ 3.30GHz 4 processor.

Figure 13 shows the MSEs of predicted energy density and second Piola-Kirchhoff stress resultants, as well as a comparison between two sampling techniques, namely uniform sampling and Sobol sampling [44], on the testing errors. Uniformly sampled data points are generated in a way that for the three features in the input vector , entries occupy the response database, where

is the number of uniformly distributed entries in each feature. In total seven datasets with

are considered; of each are taken as training points and the remaining as testing points. Moreover, is rounded to the fifth decimal number for all response databases. MSEs in stress predictions are comparatively higher than those in energy predictions but remain less than 2. Considering the error convergence and algorithmic efficiency, we use a response database with 2601 test data points in the subsequent macroscale simulations. Furthermore, this particular trained GPR model has the optimum hyperparameter values and the maximum log marginal likelihood . Model training time was recorded as 29 minutes and 48 seconds.

Figure 13: GPR model training and testing. The mean squared error of the predicted strain energy density and the second Piola-Kirchhoff stress resultants for uniform and Sobol sampling.

In Figure 16, we visualise and compare the predictions of a subsample of the chosen response database. Figure (a)a depicts the strain energy density predictions, whereas Figure (b)b presents the stress predictions for the biaxial homogenised response of an RVE.

Figure 16: GPR model training and testing. Predictions of (a) strain energy density (sorted by ascending potential) and (b) stress resultants for a biaxial strain state of , and .

5.2 Stretched membrane I: comparison of yarn-level and homogenised displacements

Figure 19: Stretched membrane I. Problem description for the yarn-level (a) and the homogenised membrane (b) models.
Figure 22: Stretched membrane I. Comparison of horizontal and vertical displacement iso-contours of the yarn-level (a) and homogenised membrane (b) models.

We consider a membrane under uniaxial tension and compare its response using the Gaussian process homogenised material model and an equivalent yarn-level model. The yarn-level model comprises 12 loops in the course direction and 12 loops in the wale direction. We use for the analysis with the homogenised material model a membrane of equal size, that is, long, wide and thick. The membrane is discretised with a structured quadrilateral mesh with elements. Problem descriptions of the two models are presented in Figure 19.

Both models are stretched by in the -direction and the resulting and displacements are presented in Figure 22. For comparison purposes, the results of the yarn level model in Figure (a)a are visualised by projecting the nodal values of the yarn onto the

plane and interpolating them on a Delaunay mesh. The

-direction displacement distribution of the yarn model in Figure (a)a is slightly different from that of the homogenised membrane in Figure (b)b. This difference has two causes. First, zero out-of-plane displacement boundary conditions on the top and bottom edges of the yarn-level model directly contribute to the observed difference. These boundary conditions are applied to simulate selvedge stitches that restrict the top and bottom yarns from being freely straightened during deformation. Secondly, the geometric asymmetry of the yarn-level model about the mid-horizontal plane has an influence on the overall response due to the relatively small number of loops used in both course and wale directions. Furthermore, a very similar Poisson effect is observed in both the yarn-level and homogenised membrane models. The absolute maximum -displacement is recorded with for a stretch by in the -direction.

5.3 Stretched membrane II: homogenised stresses

We perform two stretch tests on a square knitted membrane sheet with a side length of in the course and wale directions. A structured quadrilateral mesh of size ( elements) is used. During the displacement controlled deformation one edge is fixed while the other is stretched in the course or wale direction, respectively, by . The deformed shapes and stress contours are depicted in Figure (a)a for course-wise and in Figure (b)b for wale-wise stretching. Figure 25 clearly shows the orthotropic response of the knitted membrane as the stress resultants are different depending on the direction of the stretching, comparing, e.g.,  in Figure (a)a with  in Figure (b)b. Furthermore, the stiffer response in the wale-wise direction manifests itself in a higher stress in Figure (b)b in comparison to in Figure (a)a.

Figure 25: Stretched membrane II: Second Piola-Kirchhoff stress contours of the homogenised course-wise stretched (a) and wale-wise stretched (b) membranes.

6 Conclusions

We introduced a data-driven approach for computational homogenisation of knitted membranes to address the challenges posed by conventional homogenisation schemes. The inherent large deformations in knitted textiles are accurately captured by the finite deformation rod model and validated against experimental and numerical results. Incorporating the statistical GPR model in computational homogenisation circumvents the need for expensive microscale simulations at each quadrature point of the macroscale membrane, thus yielding significant computational savings.

The presented approach can be extended in several ways. In this paper, we considered only weft-knitted membranes with a uniform stitch pattern. To produce complex three-dimensional surfaces, it is necessary to alter the stitch pattern, for instance, by locally increasing or decreasing the number of loops from row to row or reducing the number of rows. Such changes represent discontinuities in the stitch pattern, and suitable RVEs have to be defined to capture their homogenised response. Putting aside questions of RVE size and validity of homogenisation assumptions, it is straightforward to consider in the introduced data-driven approach RVEs with other stitch patterns. Furthermore, depending on the loading conditions knitted membranes are prone to geometric instabilities in form of wrinkling on the macroscale and rod buckling on the microscale. In presence of such instabilities the homogenisation assumptions are usually not valid and alternative approaches must be used [30, 28]. The data-driven variants of these approaches are essential for the analysis of large-scale knitted membranes. Finally, a key advantage of the data-driven approach is the prospect to consider design parameters, pertaining to the stitch geometry or yarn material, in the GPR model. This opens up the possibility to optimise those parameters with an efficient gradient-based optimisation algorithm.

Appendix A Derivatives of strains and internal forces

In this Appendix we summarise the detailed equations used in the implementation of the introduced finite deformation rod finite element.

a.1 Strain derivatives

The derivative of the components  of the strain tensor (13) with respect to the nodal displacements are given by



The derivative of the rotation matrix in (9) with respect to the nodal displacements reads



The derivatives of the bending and torsional shear strains in (13) with respect to the nodal twist are given by