Review of “A simple extension of FFT-based methods to strain gradient loadings – Application to the homogenization of beams and plates with linear and non-linear behaviors”

Because of their simplicity, efficiency and ability for parallelism, FFT-based methods are very attractive in the context of numerical periodic homogenization, especially when compared to standard FE codes used in the same context. They allow applying to a unit-cell a uniform average strain with a periodic strain fluctuation that is an unknown quantity. Solving the problem allows to evaluate the complete stress-strain fields. The present work proposes to extend the use of the method from uniform loadings (i.e. uniform applied strain) to strain gradient loadings (i.e. strain fields with a uniform strain gradient) while keeping the algorithm as simple as possible. The identification of a subset of strain gradient loadings allows for a minimally invasive modification of the iterative algorithm so that the implementation is straightforward. In spite of a reduced subset of 9 independent loadings among the 18 available, the second part of the paper demonstrates that it is enough for considering the homogenization of beams and plates. A first application validates the approach and compares it to another FFT-based method dedicated to the homogenization of plates. The second application concerns the homogenization of beams, for the first time considered (to author's knowledge) with an FFT-based solver. The method applies to different beam cross-sections and the proposition of using composite voxels drastically improves the numerical solution when the beam cross-section is not conform with the spatial discretization, especially for torsion loading. As a result, the massively parallel AMITEX_FFTP code has been slightly modified and now offers a new functionality, allowing the users to prescribe torsions and flexions to beam or plate heterogeneous unit-cells.


Reviewer 1 (Anonymous)
Reviewer Note: in the following review, I will use "FFT-based solver" to denote the class of numerical method that rely on the discretization of the Lippmann-Schwinger equation (and the use of the fast Fourier transform) for the full-field simulation of heterogeneous materials in a periodic setting. This paper deals with the extension of FFT-based solvers to non-trivial prescribed macroscopic strains. The goal of the author is to propose a method which is as little intrusive as possible. This leads to considering a sub-class of strain gradients for which the iterations are indeed unchanged.
The proposed method is illustrated with a large number of test-cases, ranging from plates to beams, the latter class of problems being new within the framework of FFT-based solvers.
The paper is well organized, and the proposed method is new. It is however obscure at some places: what is prescribed and what is measured as a result is not totally clear. The differences between Ω, Ω Per and Ω * are sometimes confusing.
Authors What is prescribed and what is measured has been clarified (see answer comment 2).
Actually, the differences between Ω and Ω Per was confusing and not really relevant. It has been removed in the paper. However, the difference between Ω and Ω * is essential. Instead of solving problem (18) (now (15)) on the domain Ω, with stress free boundary conditions, problem (9) is solved on an enlarged unit-cell with full periodic boundary conditions, simulating stress free boundary conditions by additional void voxels. It has been clarified in section 3.1.1.2 and 3.2.1.2.
Reviewer Comment 1. It is not clear to me what the author is trying to achieve. The whole design of FFT solvers relies on the concept of eigenstrain or eigenstress. Hence, Problem (8) is naturally taylored to FFT-based solvers. I believe that it suffices to include * in the stress polarization as follows = ( − 0 ) :˜+ : * = − Γ 0 ( ) and the iterations should be nearly unchanged. Maybe this is still deemed too intrusive; could the author comment?
Authors I really thank the reviewer for that comment. Actually, I was focused on the definition of the polarization = : − 0 : , which led me to define additional constraints on to keep the algorithm as simple as the original. Using the definition = : − 0 :˜keeps the algorithm as simple as the original without additional constraints. I should have seen it before. . . Therefore, the code has been modified (hopefully, only slight modifications) and the description of the algorithm in section 2 has been modified (simplified!).
Reviewer Comment 2. It seems to me that there is a confusion between prescribed and measured strains.
Thus, the discussion surrounding Eqs. (18) and (19) is extremely confusing. The author mentions the Poisson effect for and which are prescribed quantities. The Poisson effect should be observed on the total strain , which solves the unit-cell problem.
Authors The description of prescribed and measured quantities has been added in Appendix A, with reference to the three problems of interest: with 3 periodicity conditions (for the enlarged unit-cells), 2 periodicity conditions (for plate unit-cell), 1 periodicity condition (for beam unit-cell).
The text between equation (18) and (19) (now (15) and (16)) has been clarified: the in-plane components of both and are prescribed, the out-of-plane components are measured. For example, for a simple test with 11 = 1, 22 = 12 = 0, the component 33 (as well as all the out-of-plane components) arising from a Poisson effect, is measured.
However, a point remains unclear to the author when applying torsion loading to the beam. According to the arguments proposed in Appendix A, the components G123 and G132 of the strain gradient are not 'prescribed'. Our numerical simulations prove that they are prescribed. Additional arguments have not been found to explain it. This point is mentioned in Appendix A3 and in section 3.2.1.
Reviewer This becomes even more confusing when the author states (without convincing arguments) that the (1, 3) and (2, 3) components of these prescribed quantities can be set to zero when adding zero-stiffness layers to the unit-cell [see around Eq. (20)].
Authors Arguments have been given in section 3.1.1.2 to explain that: • in-plane components of * and * , prescribed on Ω * , are also prescribed on Ω and, • on the other hand, the out-of-plane components of and cannot be prescribed on Ω but can be evaluated as post-treatments.
Reviewer Comment 3. I believe that problem (18) was first introduced by Caillerie and then popularized by Cecchi and Sab (2002) that ought to be cited. In this paper, the out-of-plane components (1, 3), (2, 3) and (3, 3) of the prescribed tensor are zero. This is at odds with the present paper, where the ( Authors Actually, problem (18) (now (15)) is the same as used by Cecchi and Sab (2002) as well as Nguyen, Sab, Bonnet (2008) (cited in the present paper, with an added reference to Caillerie 1984). However, where I state that the out of plane components of and on Ω should remain 'measured' quantities (see answer to comment 2), they prescribe these values to 0. The problem is that I hardly see how it can be compatible with the stress free boundary conditions · = 0 applied on the upper and lower surfaces. For example, how the average transverse strain 33 can be prescribed to 0 while, at the same time, the upper and lower surfaces are free so that 33 is free?
The following paragraph has been added.
"The author's opinion is that prescribing these quantities is not consistent with the stress free boundary condition applied on the upper and lower plate surfaces (i.e. · = 0 on 3 Ω). For example, the transverse strain 33 (but also 33 ) arising when applying a uniaxial strain 11 depends on the Poisson coefficient and can't be set to 0. As observed in section 3.2.1 the numerical results are consistent with the analytical and numerical results given in [27] for heterogeneous plates. " The problems mentioned by the reviewer in the examples have been solved (see comment below) they were coming from a technical issue, not from the proposed approach.
Reviewer Comment 4. The provided convergence diagrams are inconclusive in the present form. These plots should be presented in a log-log scale. Also, I believe that the simulations are extremely fast, so I would recommend to add more points for finer grids, for the sake of convergence test. For example, it seems that the error in Table 2 stagnates (albeit to a small value). If this is confirmed, that might indicate a bias which in turn may be caused by what I discuss in Comment 3 above.
Authors Table 2 now displays results from very low to very high spatial resolution. The new results are now conclusive: the error does not stagnate anymore. The problem was not associated to comment 3 (whose answer is given above) but to a technical problem in the input data used in the simulations. The author thanks the reviewer for this. To be short, the vtk files (3D images) used as input data are written by a matlab code and the voxel size was written in an insufficiently precise ascii format. The finer the grid size, the higher the effect of the roundoff. This problem has been corrected with a different precision of the asci format.
The results provides in Appendix have also been corrected. According to reviewers comment, the relative error plots for beam (bending and torsion) have been presented in a loglog scale (Figures 3 and 7).
Reviewer Comment 5. I believe that "flexion" should be replaced with "bending".
Authors "Flexion" has been replaced by "bending" in the paper.

Authors done
Reviewer Eq. (10), "The first constraint arises from [. . . ]". I do not think that the compatibility of * is actually required (although it can be assumed).
Authors As must be compatible (it is the solution of the problem) and as˜is compatible, I believe that * must be compatible.
Authors Yes. Added in the paper.
Reviewer Below Eq. (14): reference to index "a" is probably wrong.
Authors Yes. Moreover this paragraph has been removed.
Authors Actually, problem (6) is a part of problem (7). These expression have been replaced simply by problem (7).
Reviewer Sec. 2.5: it ought to be mentioned that the Green operator in the present case is different from the previous Green operator.

Authors Done
Reviewer Sec. 2.6: "conclusion" is weird in the middle of the manuscript. "summary"?
Authors This sub-section has been removed Reviewer Sec. 3, "his approach relied on a dedicated solver": my understanding is that this approach relies on a dedicated Green operator, the solver itself (e.g. fixed-point iterations) being unchanged.
Authors Done.
Reviewer Sec. 3.1.1, "due to the added void layers [. . . ] can be chosen arbitrarily". Why this statement is true is unclear, to me at least.
Authors See answers to comment 2 and 3 and modifications provided in the section.
Authors Important modifications have been proposed to clarify it.
Reviewer Eq. (25): the "p" superscripts are reminiscent of plastic deformations, which might be confusing in the present paper where the prescribed loading could indeed be seen as eigenstrains.
Authors Precision is given to avoid any confusion.
Reviewer The "mechanical approach" introduced in Sec. 3.1.1 should be renamed: computing the macroscopic energy is no less mechanical! Authors I understand! but I really didn't find better to distinguish the two approaches. I believe it is not so important. . .
Reviewer Sec. 3.1.3, "the decrease is less monotonous". The proposed method is supposed to be rigorous, so that the error should tend to zero. Finer grids should be considered, and convergence should be clearly analyzed in log-log space (see my previous comments).
Authors It has been removed. Decrease is now monotonous and additional results are given in Table 2 up to high spatial resolutions.
Reviewer Sec. 3.2.1, "the other components are not prescribed but do not systematically vanish". This is confusing. is a prescribed tensor, all of its components. The Poisson effect manifests itself on the total strain , I believe.
Authors This is the same comment-question than for plates. Prescribed and measured quantities have been clarified for the different problems in Appendix A. Note that no arguments are still under investigation to show that the torsion loading can be prescribed, as observed in our numerical experiments.
Authors This is the same question than for plates. The modification made for plates should clarify it.
Reviewer Between Eq. (33) and (34), "evaluated over Ω". If I understood correctly is not evaluated from the solution to the problem at hand, it is prescribed over Ω. This confusion occurs at several places in the manuscript.
Authors The point has disappeared between (33)  Authors Done Reviewer Sec. 3.2.3.1, "oscillations in the neighborhood of the beam interface". The shear stresses in fact develop a boundary layer to accommodate the boundary conditions.

Authors Done
Reviewer 2 (Sébastien Brisard) Reviewer Synopsis. This paper extends the FFT-based method of Moulinec and Suquet to strain gradient loadings and applies it to bending problems for beams and plates. This innovative and easy to implement extension has a wide range of possible applications for the characterization of composite materials. Thus, the manuscript should only be accepted provided the following MAJOR REMARKs are taken into account.
Reviewer Major remarks. I believe that the algorithm is absolutely correct (and have tested that it is working), but I have my doubts about the reasoning for the limitation of the possible strain gradient loadings: After equation (11) the author states, that the polarization is periodic. For me this is not clear. Even more, the stress should be periodic if and only if : : is periodic.
Authors These limitations have been removed since the algorithm is now written without any restriction (see answer to reviewer 1 -comment 1).
Reviewer The author should give some more details about his implementation in Section 2.4.
Authors A more detailed description of the algorithm is now provided in section 2.5.
Reviewer Section 3.1.2: The author should additionally investigate laminates with zero Poisson's ratio.
Then the assumption of the laminate theory, that the thickness of the laminate does not change would be correct.
Authors Zero Poisson's ratio is not easy to take into account as it is associated to an infinite shear modulus. A value of 10 4 has been used instead. As the main purpose of section 3.1.2 was to provide a validation from a comparison with the analytical and numerical results given in [27], these new results are provided in Appendix B.
Reviewer Section 3.1.3: The choice of constant Poisson's ratio is not realistic for fiber reinforced composites and makes the equation somewhat easier to solve. The author should investigate at least the convergence behavior for a realistic set of material parameters.
Authors Additional results are given in Appendix C, and briefly discussed in section 3.1.3, for a glass fiber/epoxy composite.
Reviewer 3. Minor remarks. Last line before Section 2.4: (with i 6=; a 6= k and k 6= i) should be replaced by (with i 6=; j 6= k and k 6= i) Authors This line was located in a part that has been removed (since the algorithm is now written without any restriction, see answer to reviewer 1 -comment 1).

Review of version 2
Permalink: hal-02942202v2 Authors Dear Editor and reviewers, Please find my second revision of the paper "A simple extension of FFT-based methods to strain gradient loadings -Application to the homogenization of beams and plates with linear and non-linear behaviors". Following Editor's advice and in order to reinforce the originality and enlarge the scope of the paper, an additional section devoted to non-linear behaviour has been added. All the points raised by the reviewers have been answered below, with corresponding modifications highlighted in blue in the paper. I sincerely thank you for your helpful remarks and hope that this modification will convince you.

Reviewer 1 (Anonymous)
Reviewer The article was substantially improved. In particular, the possible degrees of freedom for the strain gradient have been presented very nicely.
I cannot understand the problem noted in Appendix A3. Since the simulation is performed on Ω * with periodic boundary in all three space directions you can prescribe on Ω * all 9 components of the strain gradient. Since Omega and Ω * differ only by the extension with zero stiffness, all moments must act on Ω itself.
In Equation (26) is a typo: The author writes twice *

3
Since the author took all remarks of the reviewers into account the manuscript should be accepted.

Reviewer 2 (Sébastien Brisard)
Reviewer I am still not convinced of the novelty of this approach. How different is it from what Nguyen, Sab and Bonnet do in Eq. (34) and in the comments that follow that equation. Unless the author clearly asserts this point, I think that this paper should be rejected in order to be refocused around the use of zero-thickness pixel/voxel layers, which is an interesting technique (albeit not new). • Correct me if I'm wrong, but I don't think that the equivalence holds in Eq. (2). I think that if is Ω-periodic, then might be affine-plus-periodic. I believe this is discussed in the Milton book. • I agree with the author that out-of-plane components of and are not prescribed. In fact, from the point of view of plate theory, these components simply do not exist. and are in-plane tensors. • I don't see how = 0 is a problem. • All figures have extremely poor quality. • On all curves, points corresponding to actual simulations should be marked clearly. • Fig. 3 (disk+composite voxels) worries me. I don't see any convergence, here. • Fig. 5. Again, convergence should be analyzed in log-log scale. In the present case, you ought to plot log(∥ red − 1∥) as a function of log(spatial resolution). • Same goes for Fig. 9.

Editor's assessment (Julien Réthoré)
In a first review, the recommendations of the Reviewers were mitigated. Reviewer 1 was in favor of accepting the paper but Reviewer 2 was considering rejecting the paper due to a lack of novelty. In the first version of the paper, to my point of view, the novelty relied on the beam homogenization, not on the section dealing with plates. Nevertheless, compared to the paper the author is referring to for the presented theory, the proposed FFT solver is far more efficient and some of the improvements allow the author to analyse beams straightforwardly. Then in a kind of stand alone section, a strategy for finite strain was derived. This seemed to me an original part, but not very useful since only the equations were described without any numerical example. Therefore, I decided not to accept the paper directly and major revisions were asked for. In addition to a detailed response to the Reviewers' comments, the author was suggested to improve the paper, e.g. by proposing some examples for plates or beams using the proposed finite strains framework.
A second version was submitted by the author trying to enlarge the scope of the paper and to more clearly state the novelty of the proposed solution algorithm. Instead of working on large strains, the author focused on non-linear constitutive behavior. The part of the paper dealing with large strains has been removed and some applications with non-linear materials added. Further, the gap in terms of solution algorithm compared to the reference by Nguyen et al. is now more clearly established. Even if the theoretical aspects are similar to what was proposed by Nguyen et al., the ability of the strategy described in the present paper to deal with arbitrary cross-section is an important point. While Reviewer 2 recommends to reject the paper, I think the author now better emphasizes on the novelty of his strategy which increases significantly the scope of applications of the underlying theoretical formulation. Further, the validation of the proposed implementation assessed by reproducing most of the original results by Nguyen et al. provides a convincing demonstration of its robustness. The revised version also contained clear answers to all comments of the reviewers. I thus decided to accept the revised version of the paper.