Overall Thoughts

Overall, this is a well-formulated paper with clear objectives. The authors propose a method where the boundary constraints in a layout optimization problem are implemented as inequalities for the distance of the turbine to the nearest edge of the boundary polygon. To avoid the issue with these distinct inclusion zones—namely that the initial conditions fix the number of turbines in each zone—the authors propose two methods: (1) relaxing the wind farm boundaries for a certain number of iterations at the beginning of the optimization and (2) greedily placing the turbines in feasible locations (where the AEP of the new turbine is maximized) before starting the optimization. My major comments are about the explanation of the relaxation process in the “Approach 2” and about the comparison of this method to other penalty methods in the literature. Other minor comments about the text are also listed below. I would recommend this paper for publication in Wind Energy Science with these minor revisions.

Major Comments

Penalty methods: The introduction describes how penalty methods have been used to handle boundary constraints in other papers, but the downside is the tunable hyperparameter. Approach 1 for the optimization uses no penalty methods, and the major downside is that the turbines are fixed in their respective inclusion zones after their initial placement. If a penalty parameter were included in Approach 1, wouldn’t the turbines be able to jump between the inclusion zones more easily (as discussed around Line 402), and the algorithm would perform better? I suppose that the convergence speed would suffer as a result, but I would expect the solutions to compare better to Approaches 2 and 3.

Boundary relaxation: I find Section 3.4 confusing to digest. Starting at Line 317, what is the importance of the “last 300 m?” Does 300 m force the inclusion zones to combine with one another, similar to the first panel in Fig. 3? Second, does the un-relaxation at the rate k_r start at 100 iterations and continue until the limit \gamma_r is reached? If so, how does Eq. 17 reconcile with Eq. 9? Or is the return to the actual boundaries being done continuously at the rate k_r starting at the first iteration? Lastly, on Line 324: if k_r = 100, \gamma_r should equal 103, not 97?

Line 298 mentions how “The values of k_r and \gamma_r have to be selected accordingly to the size of the domain and the inclusion zone areas.” This sentence is vague to me: for a larger domain with more turbines, would the number of iterations increase? If the inclusion zone areas are more complex, would these values increase? Are there heuristics for how many iterations the optimization will take such that \gamma_r is chosen to be smaller than that expected limit? Also, I think calling k_r a relaxation speed (Line 319) is confusing because it’s really an un-relaxation speed. There’s surely a better name than that but ‘relaxation’ implies that the relaxation is progressing with k_r, not being undone (i.e. boundaries are returning from their relaxed state to the normal state).

Lastly, Line 347 discusses the differences between the performance of Combinations 3 and 4. The authors explain this by saying that Combination 4 does not have enough “time to find good positions”, but I think this really has to do with the size of the relaxed inclusion areas. The number of relaxation iterations multiplied by the relaxation speed is the maximum offset of the boundary, and with the smaller k_r and fixed \gamma_r, the relaxed boundaries are not sufficiently large to really explore the unconstrained solution space.

Minor Comments

Line 74: Can more detail be given on why gradient-based optimization would be ineffective? I suppose you’re inferring that the high sensitivity to the initial conditions would require the user to run many optimization studies with randomized initial conditions to find the best candidates, which would be more inefficient than the proposed methods.

Line 78: “This is achieved by adding an offset to the distance determined by the method.” This sentence is vague—adding an offset to the distance to the boundary of the polygon?

Eq 1: Maybe a preface about the lower and upper limits of the design variables could be included here. In Fig. 3, it is clear why these limits on the computational domain are necessary, but from Eq. 1 it appears that the function C_i makes them redundant (since the reader assumes that the inclusion zone boundaries are within the computational domain).

Eq. 4: Can more detail be given on why this step is necessary, as opposed to just using the sign of the distance to the nearest vertex (a_{ik} or b_{ik})? It seems that you’d just want to move towards that nearest point on the boundary. Also, where do the ‘k’ and ‘k-1’ indices come from? Everything else in this section is with respect to indices ‘k’ and ‘k+1’.

Eq. 6: A sign is needed here to denote which side of the boundary the turbine is on. Right now, the absolute value makes C_i >= 0 always, which would satisfy the constraint even when the turbine is outside of the inclusion zone.

Eq. 11: The subscript ‘u’ should be defined similar to the subscript ‘d’. Also, are the authors assuming that the speedup factor does not depend on the wind speed? That seems reasonable, but the assumption could be stated explicitly.

Line 254: Why does the crosswind distance not follow the terrain like the upstream distance does?

Line 270: AEP is from Section 2.4, not Equation 2.4, and the wind speed and direction are defined in Eqs. 11 and 12, not 12 and 13.

Line 275: 12 wind sectors are mentioned here, but then Line 284 mentions 1 degree wind direction resolution for the optimization. What is the meaning of the 12 wind sectors?

Line 290: Can this “chunkification” be explained more? Does this mean that the wind direction resolution is not 1 degree, or is this something to do with how the vectorized wind conditions are distributed across the computational nodes?

Line 294: Autograd is used for automatic differentiation of the objective function and also the spacing constraint? I assume that the analytic gradients of the boundary constraint function (Eq. 7 and 8) have been provided to the optimizer.

Line 326: I would be careful to clarify between optimization “time” and “iterations” throughout the paper. At some points, “time” is used where the authors are really referring to iterations and convergence.

Line 329: Is there any explanation why these seeds fail beyond that the solutions violate the constraints? I suppose it’s tied to the fact that the longer relaxation periods with slower relaxation speeds allow the optimizer to settle into local optima that are infeasible, but maybe that could be explained a bit more again.

Line 364: What does “incompatibility of inequality constraints” mean?

Line 378: The lower cost of the smart-start optimizations, even including the time for the initialization (because the convergence is faster), is a neat result. Is there a sensitivity of the smart-start results to the choice of the randomness parameter? How much randomness is required when initially placing the turbines to compensate for the unknown potential optimal AEP in each inclusion zone based on that number of turbines?

Line 415: What is meant by “vectorized” here? Are the wind conditions computed in PyWake in a vectorized fashion? Are the constraints implemented in TOPFARM in a vectorized fashion? Vectorization has not been mentioned up until this point.

Technical Comments

Line 91: Small typo: “TOPFARMS relies [on] PyWake”

Line 243: Symbol for U_\infty is different from the other notation.

Line 341: Typo, remove “for [a] the”.

The manuscript in question describes layout optimization of an offshore wind farm based on maximizing AEP (using a simple Gaussian wake model), under consideration of inclusion or exclusion zones. Instead of commonly employed heuristic approaches state of the art gradient-based optimization is used. For this it is necessary to determine gradients of the constraints. The authors solve this problem by calculating the (signed) distance of each turbine to the nearest inclusion / exclusion zone and taking its derivative. A second issue that is addressed in the paper is the problem of many local extrema. In particular, a turbine that starts in an inclusion zone is never able to move to another inclusion zone, and the number of turbines per inclusion zone is a constant. To address this problem, the authors try two approaches. In the first approach, the distance constraint is initially relaxed and becomes closer to the true constraint with each iteration. In the second approach, called "smart start", an initial positioning of the turbines is calculated based on the simple heuristic of placing wind turbines sequentially at the remaining most promising location.  This seems to be a so-called greedy algorithm.

Results of a case study are presented that show a surprisingly large variability in "optimal" configurations after convergence, and that while the more classical relaxation approach improves the solution on the average, the smart start strategy performs even better.

Overall the paper is concisely and well written. The idea to use a distance function, although relatively straightforward, is carefully developed and explained. Although results are somewhat worse than I would have expected, the technical developments and their discussion seems to me to be valuable to readers of Wind Energy Science, and I am quite in favor of publishing this manuscript.

That said, I have a number of comments on the paper that I think could make the paper even more interesting to readers, so I would like the authors to comment on them. 

Detailed comments:

1. line 50ff: Why are approaches using polygons or penalty functions not sufficient, in the opinion of the authors? A few details would be welcome here.
2. line 118f: "If more than one polygon of the same type overlap, these are merged into a single one". This is fine, in principle. But what about inclusion and exclusion zones - I assume that these are not allowed to overlap?
3. line 136: What is a "unitary" vector? I assume the authors mean that these vectors are normalized to unit length?
4. Figs. 1-2: It is very confusing to see the normal vectors drawn at the starting node of their respective edge, instead of - as I would expect - somewhere around the middle of the edge. Mathematically, there is no problem, but it goes somehow against expectations of how to visualize these things. Please reconsider.
5. line 218: "until all wind turbines are added." - I would assume that this "greedy" algorithm could also fail under certain conditions? Moreover, since it seems to be deterministic (when r = 0), there would then be no starting condition for the actual optimization? Is this the (main) reason for adding randomness to this procedure? (I think there are straightforward ways to deal with this problem better, though - maybe something the authors might want to think about)
6. Eq. 15: I assume the gradients of the velocity deficit are also needed, when determining the gradient of the AEP? But this is handled by algorithmic differentiation, correct?
7. line 275: Typo: "Weibull ditribution"
8. line 306: "The idea ... is to ... prove that the nearest distance method succeeds..." - Well, "proving" this might be asking for too much. Why not say "is to ... demonstrate"?
9. line 314ff: "there might not be enough time (iterations) to explore the domain". Well, the relaxation of the distance by way of Eq. 9 prescribes a specific way in which this constraint is changing from iteration to iteration. In some cases the optimization even terminates when the distance function is still relaxed (line 361), which the authors consider to be invalid optimization runs!? The first seems somewhat ad-hoc, and the second seems somewhat strange. I would imagine that it would be much more natural to fully optimize the wind farm layout for a fixed offset of the distance, and only then adjust the constraint and solve the next (slightly less relaxed) optimization problem - using the previous optimal solution as an initial (though most likely infeasible) state. Thus, instead of considering one optimization problem with continuously changing constraint relaxation this would be a series of optimization problems with different (but fixed) constraint relaxations - until the final problem with no relaxation is solved. Why did the authors not try this, very often used, strategy?
10. Continuing from the previous comment, there might also be better ways to change the relaxation offset to the distance adaptively, i.e., depending on the speed of convergence (how much better the solution becomes between iterations).
11. line 321: "before relaxation is applied to the last 300 meters." - There seems to be some confusion about the use of the term "relaxation" and I would like the authors to review their use of it. The relaxed optimization problem has the least strict constraint. It is thus not correct to talk about relaxing the constraint in later iterations where the constraint actually becomes more strict! Similar issues of imprecise use occur in line 346 ("no clear benefit in slowing down the relaxation") and line 375 ("high speed of the relaxation").
12. line 325: Why the value of 97? Why not 103? There seems to be some inconsistency between this here and Eq. 17.
13. Table 2: A comment seems appropriate why not both relaxation and smart-start have been tested together.
14. page 17. The results are discussed mainly in terms of average AEP, but what about the maximum values? These are, after all, what one is looking for when trying to optimize the problem!
15. Fig. 9: I am somewhat surprised to see so much variability in the results, depending on the used random seed / initial position of the turbines. Obviously there seem to be many local maxima in this problem. Although the approach of the authors is promising, this seems still to be a major problem - how to find the actual, global optimum? Some thoughts along this direction would be highly appreciated in the discussion.
16. Coming back to the problem of local maxima, have the authors considered to use the idea of marginal improvements - i.e., estimating numerically how much the addition or removal of a single wind turbine from each inclusion zone would add / remove from the overall AEP (per turbine)? This type of sensitivity could then be used to determine when or how to move turbines between inclusion zones...