# Yamabe systems and optimal partitions on manifolds with symmetries

Mónica Clapp111M. Clapp was partially supported by CONACYT (Mexico) through
the grant for project A1-S-10457. and Angela Pistoia222A. Pistoia was
partially supported by Fondi di Ateneo “Sapienza” Università di Roma (Italy).

###### Abstract

We prove the existence of regular optimal $G$-invariant partitions, with an
arbitrary number $\ell\geq 2$ of components, for the Yamabe equation on a
closed Riemannian manifold $(M,g)$ when $G$ is a compact group of isometries
of $M$ with infinite orbits. To this aim, we study a weakly coupled
competitive elliptic system of $\ell$ equations, related to the Yamabe
equation. We show that this system has a least energy $G$-invariant solution
with nontrivial components and we show that the limit profiles of the its
components separate spatially as the competition parameter goes to $-\infty$,
giving rise to an optimal partition. For $\ell=2$ the optimal partition
obtained yields a least energy sign-changing $G$-invariant solution to the
Yamabe equation with precisely two nodal domains.

Keywords: Competitive elliptic system, Riemannian manifold, critical
nonlinearity, optimal partition, free boundary problem, regularity, Yamabe
equation, sign-changing solution.

MSC2020: 35B38, 35J20, 35J47, 35J60, 35R35, 49K20, 49Q10, 58J05.

Dedicated to Norman Dancer on the occasion of his 75th birthday

##  1 Introduction

This paper is concerned with the existence and asymptotic behavior of
solutions to the weakly coupled competitive Yamabe system

$\begin{cases}\mathscr{L}_{g}u_{i}:=-\Delta_{g}u_{i}+\kappa_{m}S_{g}u_{i}=|u_{i}|^{2^{*}-2}u_{i}+\sum\limits_{\begin{subarray}{c}j=1\\\
j\neq
i\end{subarray}}^{\ell}\lambda_{ij}|u_{j}|^{\frac{2^{*}}{2}}|u_{i}|^{\frac{2^{*}}{2}}u_{i}\quad\text{on
}M,\\\ u_{i}\neq 0,\qquad i=1,\ldots,\ell,\end{cases}$ (1.1)

where $(M,g)$ is a closed Riemannian manifold of dimension $m\geq 3$, $S_{g}$
is its scalar curvature, $\Delta_{g}:=\mathrm{div}_{g}\nabla_{g}$ is the
Laplace-Beltrami operator, $\kappa_{m}:=\frac{m-2}{4(m-1)}$,
$2^{*}:=\frac{2m}{m-2}$ is the critical Sobolev exponent, and
$\lambda_{ij}=\lambda_{ji}<0$. We assume that the quadratic form induced by
the conformal Laplacian $\mathscr{L}_{g}$ is coercive.

This system was recently studied by Clapp, Pistoia and Tavares [9]. Here we
complement the results obtained in [9] by considering manifolds with
symmetries.

Let $G$ be a compact group of isometries of $M$ and let $Gp:=\\{\gamma
p:\gamma\in G\\}$ denote the $G$-orbit of a point $p\in M$. Recall that a
subset $X$ of $M$ is said to be $G$-invariant if $Gp\subset X$ for every
$p\in\Omega$ and a function $u:X\to\mathbb{R}$ is called $G$-invariant if it
is constant on every $G$-orbit of $\Omega$. We shall say that a solution
$(u_{1},\ldots,u_{\ell})$ to the system (1.1) is $G$-invariant if every
component $u_{i}$ is $G$-invariant.

We prove the following result.

######  Theorem 1.1.

If $1\leq\dim(Gp)<m$ for every $p\in M$, then the system (1.1) has a least
energy $G$-invariant solution and infinitely many $G$-invariant solutions.

To describe the limit profile of least energy $G$-invariant solutions to the
system (1.1) as $\lambda_{ij}\to-\infty$, we consider the Dirichlet problem

$\begin{cases}\mathscr{L}_{g}u=|u|^{2^{*}-2}u&\text{ in }\Omega,\\\ u=0&\text{
on }\partial\Omega,\\\ u\text{ is }G\text{-invariant},\end{cases}$ (1.2)

in a $G$-invariant open subset $\Omega$ of $M$, and set

$c_{\Omega}^{G}:=\inf\Big{\\{}\frac{1}{m}\int_{M}|u|^{2^{*}}\mathrm{\,d}\mu_{g}:u\neq
0,\ u\text{ solves }\eqref{eq:u}\Big{\\}}.$ (1.3)

Let
$\mathcal{P}_{\ell}^{G}:=\\{\\{\Omega_{1},\ldots,\Omega_{\ell}\\}:\,\Omega_{i}\neq\emptyset\text{
is a }G\text{-invariant open subset of }M\text{ and
}\Omega_{i}\cap\Omega_{j}=\emptyset\text{ if }i\neq j\\}$. We shall say that
$\\{\Omega_{1},\ldots,\Omega_{\ell}\\}\in\mathcal{P}_{\ell}^{G}$ is _an
optimal $(G,\ell)$-partition_ for the Yamabe equation

$\mathscr{L}_{g}u=|u|^{2^{*}-2}u\qquad\text{on }M,$ (1.4)

if

$\sum_{i=1}^{\ell}c_{\Omega_{i}}^{G}=\inf_{\\{\Theta_{1},\ldots,\Theta_{\ell}\\}\in\mathcal{P}_{\ell}^{G}}\;\sum_{i=1}^{\ell}c_{\Theta_{i}}^{G}.$

The relation between variational elliptic systems having large competitive
interaction and optimal partition problems was observed by Conti, Terracini
and Verzini in [12, 13] and Chang, Lin, Lin and Lin in [5], and has been
extensively studied since. An ample list of references is given in [9]. Our
next result describes this relation for the system (1.1).

######  Theorem 1.2.

Assume $1\leq\dim(Gp)<m$ for every $p\in M$. Let $\lambda_{n}<0$ be such that
$\lambda_{n}\to-\infty$. For each $n\in\mathbb{N}$, let
$(u_{n,1},\ldots,u_{n,\ell})$ be a least energy $G$-invariant solution to the
system (1.1) with $\lambda_{ij}=\lambda_{n}$ for all $i\neq j$, such that
$u_{n,i}>0$ for all $n\in\mathbb{N}$ and $i=1,\ldots\ell$. Then, after passing
to a subsequence, we have that

  * $(i)$

$u_{n,i}\to u_{\infty,i}$ strongly in
$H_{g}^{1}(M)\cap\mathcal{C}^{0,\alpha}(M)$ for every $\alpha\in(0,1)$, where
$u_{\infty,i}\geq 0$, $u_{\infty,i}\neq 0$, and $u_{\infty,i}|_{\Omega_{i}}$
is a least energy $G$-invariant solution to the problem (1.2) in
$\Omega_{i}:=\\{p\in M:u_{\infty,i}(p)>0\\}$ for each $i=1,\ldots,\ell$.
Moreover,

$\int_{M}\lambda_{n}u_{n,i}^{\frac{2^{*}}{2}}u_{n,j}^{\frac{2^{*}}{2}}\to
0\text{ \ as \ }n\to\infty\quad\text{whenever }i\neq j.$

  * $(ii)$

$u_{\infty,i}$ is $G$-invariant and $u_{\infty,i}\in\mathcal{C}^{0,1}(M)$ for
each $i=1,\ldots,\ell$.

  * $(iii)$

$\\{\Omega_{1},\ldots,\Omega_{\ell}\\}\in\mathcal{P}_{\ell}^{G}$ and it is an
optimal $(G,\ell)$-partition for the Yamabe equation on $(M,g)$. In
particular, each $\Omega_{i}$ is connected.

  * $(iv)$

$M\smallsetminus\bigcup_{i=1}^{\ell}\Omega_{i}=\mathscr{R}\cup\mathscr{S}$,
where $\mathscr{R}\cap\mathscr{S}=\emptyset$, $\mathscr{R}$ is an
$(m-1)$-dimensional $\mathcal{C}^{1,\alpha}$-submanifold of $M$ and
$\mathscr{S}$ is a closed subset of $M$ with Hausdorff measure $\leq m-2$. In
particular, $M=\cup_{i=1}^{\ell}\overline{\Omega}_{i}$. Moreover, given
$p_{0}\in\mathscr{R}$ there exist $i\neq j$ such that

$\lim_{p\to p_{0}^{+}}|\nabla_{g}u_{i}(p)|^{2}=\lim_{p\to
p_{0}^{-}}|\nabla_{g}u_{j}(p)|^{2}\neq 0,$

where $p\to p_{0}^{\pm}$ are the limits taken from opposite sides of
$\mathscr{R}$, and for $p_{0}\in\mathscr{S}$ we have

$\lim_{p\to p_{0}}|\nabla_{g}u_{i}(p)|^{2}=0\quad\text{for every
}i=1,\ldots,\ell.$

  * $(v)$

If $\ell=2$, then $u_{\infty,1}-u_{\infty,2}$ is a least energy $G$-invariant
sign-changing solution to the Yamabe equation (1.4).

In the nonsymmetric case, the existence of a solution to the system (1.1) was
recently shown in [9, Theorem 1.2] under the additional assumptions that
$(M,g)$ has either dimension $3$ and is not conformal to the standard
$3$-sphere, or $(M,g)$ is not locally conformally flat and $\dim M\geq 9$. A
nonsymmetric version of Theorem 1.2 was also established in that work,
provided $(M,g)$ is not locally conformally flat, $\dim M\geq 10$, and
satisfies a further geometric condition whenever $\dim M=10$; see [9, Theorem
1.2].

The existence of a least energy sign-changing solution to the Yamabe equation
(1.4) was established by Ammann and Humbert in [1] for manifolds of dimension
$\geq 11$ that are not locally conformally flat, without any symmetry
assumption, whereas existence of infinitely many $G$-invariant sign-changing
solutions to this equation was shown by Clapp and Fernández [6] for $G$ as
above.

An immediate consequence of Theorems 1.1 and 1.2 is the following result.

######  Theorem 1.3.

If $1\leq\dim(Gp)<m$ for every $p\in M$, then for every $\ell\geq 2$ there
exists an optimal $(G,\ell)$-partition $\\{\Omega_{1},\ldots,\Omega_{\ell}\\}$
for the Yamabe equation on $(M,g)$ with the properties stated in items $(iii)$
and $(iv)$ of Theorem 1.2.

When $M$ is the standard sphere $\mathbb{S}^{m}$ and $G$ is suitably chosen,
the optimal partition admits a more accurate description as shown by Clapp and
Pistoia [7] and Clapp, Saldaña and Szulkin [10] . We shall give an account of
known results for the standard sphere in Section 4 of this paper. Section 2 is
devoted to the proof of Theorem 1.1 and Section 3 to that of Theorem 1.2.

##  2 Existence of multiple solutions

Let $H^{1}_{g}(M)$ be the Sobolev space of square integrable functions on $M$
having square integrable first weak derivatives. We write
$\langle\,\cdot\,,\,\cdot\,\rangle$ for the Riemannian metric in $(M,g)$ and
denote the induced norm by $|\,\cdot\,|$. For $u,v\in H^{1}_{g}(M)$ let

$\langle
u,v\rangle_{g}:=\int_{M}\left(\langle\nabla_{g}u,\nabla_{g}v\rangle+\kappa_{m}S_{g}uv\right)\mathrm{\,d}\mu_{g}\qquad\text{and}\qquad\|u\|_{g}:=\sqrt{\langle
u,u\rangle_{g}},$

where $\nabla_{g}$ is the weak gradient. Since we are assuming that the
conformal Laplacian $\mathscr{L}_{g}$ is coercive, $\|\cdot\|_{g}$ is a norm
in $H^{1}_{g}(M)$, equivalent to the standard one.

Let $G$ be a closed subgroup of the group of isometries of $M$ and set

$H^{1}_{g}(M)^{G}:=\\{u\in H^{1}_{g}(M):u\text{ is }G\text{-invariant}\\}.$

If $\dim(Gp)<m$ for $p\in M$ then, for any given $k\in\mathbb{N}$, we may
choose $u_{1},\ldots u_{k}\in H^{1}_{g}(M)^{G}$ having pairwise disjoint
supports. Hence, $H^{1}_{g}(M)^{G}$ has infinite dimension. The following
result will play a crucial role in the proof of Theorem 1.1.

######  Lemma 2.1.

If $\dim(Gp)\geq 1$ for every $p\in M$, then the embedding
$H^{1}_{g}(M)^{G}\hookrightarrow L_{g}^{2^{*}}(M)$ is compact.

###### Proof.

See [24, Corollary 1]. ∎

We assume from now on that $1\leq\dim(Gp)<m$ for every $p\in M$. Set
$\mathcal{H}^{G}:=(H^{1}_{g}(M)^{G})^{\ell}$ and let
$\mathcal{J}:\mathcal{H}^{G}\to\mathbb{R}$ be given by

$\displaystyle\mathcal{J}(u_{1},\ldots,u_{\ell}):=$
$\displaystyle\frac{1}{2}\sum_{i=1}^{\ell}\|u_{i}\|_{g}^{2}-\frac{1}{2^{*}}\sum_{i=1}^{\ell}\int_{M}|u_{i}|^{2^{*}}\mathrm{\,d}\mu_{g}-\frac{1}{2^{*}}\mathop{\sum_{i,j=1}^{\ell}}_{j\neq
i}\int_{M}\lambda_{ij}|u_{j}|^{\frac{2^{*}}{2}}|u_{i}|^{\frac{2^{*}}{2}}\mathrm{\,d}\mu_{g}.$
(2.1)

This functional is of class $\mathcal{C}^{1}$ and its partial derivatives are

$\displaystyle\partial_{i}\mathcal{J}(u_{1},\ldots,u_{\ell})v=$
$\displaystyle\,\langle
u_{i},v\rangle_{g}-\int_{M}|u_{i}|^{2^{*}-2}u_{i}v\mathrm{\,d}\mu_{g}-\mathop{\sum_{j=1}^{\ell}}_{j\neq
i}\int_{M}\lambda_{ij}|u_{j}|^{\frac{2^{*}}{2}}|u_{i}|^{\frac{2^{*}}{2}-2}u_{i}v\mathrm{\,d}\mu_{g}.$

So, by the principle of symmetric criticality [26], the critical points of
$\mathcal{J}:\mathcal{H}^{G}\to\mathbb{R}$ are the $G$-invariant solutions to
the system

$\mathscr{L}_{g}u_{i}=|u_{i}|^{2^{*}-2}u_{i}+\sum\limits_{\begin{subarray}{c}j=1\\\
j\neq
i\end{subarray}}^{\ell}\lambda_{ij}|u_{j}|^{\frac{2^{*}}{2}}|u_{i}|^{\frac{2^{*}}{2}}u_{i}\quad\text{on
}M,\qquad i=1,\ldots,\ell.$

We are interested in solutions $(u_{1},\ldots,u_{\ell})$ such that every
$u_{i}$ is nontrivial. They belong to the set

$\mathcal{N}^{G}:=\\{(u_{1},\ldots,u_{\ell})\in\mathcal{H}^{G}:u_{i}\neq
0,\;\partial_{i}\mathcal{J}(u_{1},\ldots,u_{\ell})u_{i}=0,\;\forall
i=1,\ldots,\ell\\}.$ (2.2)

Note that
$\mathcal{J}(u_{1},\ldots,u_{\ell})=\frac{1}{m}\sum_{i=1}^{\ell}\|u_{i}\|_{g}^{2}$
if $(u_{1},\ldots,u_{\ell})\in\mathcal{N}^{G}$.

######  Lemma 2.2.

There exists $d_{0}>0$ such that $\|u_{i}\|_{g}^{2}\geq d_{0}$ for every
$(u_{1},\ldots,u_{\ell})\in\mathcal{N}^{G}$ and $i=1,\ldots,\ell$. Therefore,
$\mathcal{N}^{G}$ is a closed subset of $\mathcal{H}^{G}$ and
$\inf_{\mathcal{N}^{G}}\mathcal{J}>0$.

###### Proof.

As $\lambda_{ij}<0$, the Sobolev embedding yields a positive constant $C$ such
that

$\displaystyle\|u_{i}\|_{g}^{2}\leq\int_{M}|u_{i}|^{2^{*}}\mathrm{\,d}\mu_{g}\leq
C\|u_{i}\|_{g}^{2^{*}}\quad\text{ for every \
}(u_{1},\ldots,u_{\ell})\in\mathcal{N}^{G},\ i=1,\ldots,\ell,$

and our claims follow. ∎

A solution $\bar{u}\in\mathcal{N}^{G}$ to the system (1.1) satisfying
$\mathcal{J}(\bar{u})=\inf_{\mathcal{N}^{G}}\mathcal{J}$ is called a _least
energy $G$-invariant solution_.

The variational approach introduced in [11] can be immediately adapted to
establish the existence of infinitely many fully nontrivial critical points of
$\mathcal{J}$. We sketch this procedure.

Given $\bar{u}=(u_{1},\ldots,u_{\ell})$ and
$\bar{s}=(s_{1},\ldots,s_{\ell})\in(0,\infty)^{\ell}$, we write
$\bar{s}\bar{u}:=(s_{1}u_{1},\ldots,s_{\ell}u_{\ell})$. Let
$\mathcal{S}^{G}:=\\{u\in H_{g}^{1}(M)^{G}:\|u\|_{g}=1\\}$,
$\mathcal{T}^{G}:=(\mathcal{S}^{G})^{\ell}$, and

$\mathcal{U}^{G}:=\\{\bar{u}\in\mathcal{T}^{G}:\bar{s}\bar{u}\in\mathcal{N}^{G}\text{
\ for some \ }\bar{s}\in(0,\infty)^{\ell}\\}.$

######  Lemma 2.3.

  * $(i)$

Let $\bar{u}\in\mathcal{T}^{G}$. If there exists
$\bar{s}_{\bar{u}}\in(0,\infty)^{\ell}$ such that
$\bar{s}_{\bar{u}}\bar{u}\in\mathcal{N}^{G}$, then $\bar{s}_{\bar{u}}$ is
unique and satisfies

$\mathcal{J}(\bar{s}_{\bar{u}}\bar{u})=\max_{\bar{s}\in(0,\infty)^{\ell}}\mathcal{J}(\bar{s}\bar{u}).$

  * $(ii)$

$\mathcal{U}^{G}$ is a nonempty open subset of the Hilbert manifold
$\mathcal{T}^{G}$.

  * $(iii)$

The map $\mathcal{U}^{G}\to\mathcal{N}^{G}$ given by
$\bar{u}\mapsto\bar{s}_{\bar{u}}\bar{u}$ is a homeomorphism.

###### Proof.

The proof is exactly the same as that of [11, Proposition 3.1]. ∎

Define $\Psi:\mathcal{U}^{G}\to\mathbb{R}$ by

$\Psi(\bar{u}):=\mathcal{J}(\bar{s}_{\bar{u}}\bar{u}).$

If $\Psi$ is of class $\mathcal{C}^{1}$ the norm of $\Psi^{\prime}(\bar{u})$
in the cotangent space $\mathrm{T}_{\bar{u}}^{*}(\mathcal{T}^{G})$ to
$\mathcal{T}^{G}$ at $\bar{u}$ is defined as

$\|\Psi^{\prime}(\bar{u})\|_{*}:=\sup\limits_{\begin{subarray}{c}\bar{v}\in\mathrm{T}_{\bar{u}}(\mathcal{U}^{G})\\\
\bar{v}\neq
0\end{subarray}}\frac{|\Psi^{\prime}(\bar{u})\bar{v}|}{\|\bar{v}\|_{g}},$

where $\mathrm{T}_{\bar{u}}(\mathcal{U}^{G})$ is the tangent space to
$\mathcal{U}^{G}$ at $\bar{u}$.

A sequence $(\bar{u}_{n})$ in $\mathcal{U}^{G}$ is called a
$(PS)_{c}^{G}$_-sequence for_ $\Psi$ if $\Psi(\bar{u}_{n})\to c$ and
$\|\Psi^{\prime}(\bar{u}_{n})\|_{*}\to 0$, and $\Psi$ is said to satisfy the
$(PS)_{c}^{G}$_-condition_ if every such sequence has a convergent
subsequence. Similarly, a $(PS)_{c}^{G}$_-sequence for_ $\mathcal{J}$ is a
sequence $(\bar{u}_{n})$ in $\mathcal{H}^{G}$ such that
$\mathcal{J}(\bar{u}_{n})\to 0$ and
$\|\mathcal{J}^{\prime}(\bar{u}_{n})\|_{(\mathcal{H}^{G})^{\prime}}\to 0$.
$\mathcal{J}$ satisfies the $(PS)_{c}^{G}$_-condition_ if any such sequence
has a convergent subsequence. As usual, $(\mathcal{H}^{G})^{\prime}$ stands
for the dual space of $\mathcal{H}^{G}$.

######  Lemma 2.4.

  * $(i)$

$\Psi\in\mathcal{C}^{1}(\mathcal{U}^{G},\mathbb{R})$,

$\Psi^{\prime}(\bar{u})\bar{v}=\mathcal{J}^{\prime}(\bar{s}_{\bar{u}}\bar{u})[\bar{s}_{\bar{u}}\bar{v}]\quad\text{for
all }\bar{u}\in\mathcal{U}^{G}\text{ and
}\bar{v}\in\mathrm{T}_{\bar{u}}(\mathcal{U}^{G}),$

and there exists $d_{0}>0$ such that

$d_{0}\,\|\mathcal{J}^{\prime}(\bar{s}_{\bar{u}}\bar{u})\|_{(\mathcal{H}^{G})^{\prime}}\leq\|\Psi^{\prime}(\bar{u})\|_{*}\leq|\bar{s}_{\bar{u}}|_{\infty}\|\mathcal{J}^{\prime}(\bar{s}_{\bar{u}}\bar{u})\|_{(\mathcal{H}^{G})^{\prime}}\quad\text{for
all }\bar{u}\in\mathcal{U}^{G},$

where $|\bar{s}|_{\infty}=\max\\{|s_{1}|,\ldots,|s_{q}|\\}$ if
$\bar{s}=(s_{1},\ldots,s_{q})$.

  * $(ii)$

Let $\bar{u}_{n}\in\mathcal{U}^{G}$. If $(\bar{u}_{n})$ is a
$(PS)_{c}^{G}$-sequence for $\Psi$, then $(\bar{s}_{\bar{u}_{n}}\bar{u}_{n})$
is a $(PS)_{c}^{G}$-sequence for $\mathcal{J}$.

  * $(iii)$

Let $\bar{u}\in\mathcal{U}^{G}$. Then, $\bar{u}$ is a critical point of $\Psi$
if and only if $\bar{s}_{\bar{u}}\bar{u}$ is a critical point of $\mathcal{J}$
if and only if $\bar{s}_{\bar{u}}\bar{u}$ is a $G$-invariant solution of
(1.1).

  * $(iv)$

If $(\bar{u}_{n})$ is a sequence in $\mathcal{U}^{G}$ and
$\bar{u}_{n}\to\bar{u}\in\partial(\mathcal{U}^{G})$, then
$\Psi(\bar{u}_{n})\to\infty$.

  * $(v)$

$\bar{u}\in\mathcal{U}^{G}$ if and only if $-\bar{u}\in\mathcal{U}^{G}$, and
$\Psi(\bar{u})=\Psi(-\bar{u})$.

###### Proof.

The proof of these statements is exactly the same as that of [11, Theorem
3.3]. ∎

######  Lemma 2.5.

$\Psi$ satisfies the $(PS)_{c}^{G}$-condition for every $c\in\mathbb{R}$.

###### Proof.

Let $(\bar{v}_{n})$ be a $(PS)_{c}^{G}$-sequence for $\mathcal{J}$ with
$\bar{v}_{n}\in\mathcal{N}^{G}$. A standard argument shows that this sequence
is bounded. Then, using Lemma 2.1 as in [7, Proposition 3.6], one sees that
$(\bar{v}_{n})$ contains a convergent subsequence. The claim now follows from
Lemmas 2.4$(ii)$ and 2.3$(iii)$. ∎

Let $\mathcal{Z}$ be a nonempty subset of $\mathcal{T}^{G}$ such that
$\bar{u}\in\mathcal{Z}$ if and only if $-\bar{u}\in\mathcal{Z}$. Recall that
the _genus of $\mathcal{Z}$_, denoted $\mathrm{genus}(\mathcal{Z})$, is the
smallest integer $k\geq 1$ such that there exists an odd continuous function
$\mathcal{Z}\rightarrow\mathbb{S}^{k-1}$ into the unit sphere
$\mathbb{S}^{k-1}$ in $\mathbb{R}^{k}$. If no such $k$ exists, we define
$\mathrm{genus}(\mathcal{Z}):=\infty$. We set $\mathrm{genus}(\emptyset):=0$.

######  Lemma 2.6.

$\mathrm{genus}(\mathcal{U}^{G})=\infty$.

###### Proof.

This is shown following the argument in [11, Lemma 4.5]. ∎

###### Proof of Theorem 1.1.

Lemma 2.4$(iv)$ implies that $\mathcal{U}^{G}$ is positively invariant under
the negative pseudogradient flow of $\Psi$, so the usual deformation lemma
holds true for $\Psi$ in $\mathcal{U}^{G}$, see e.g. [27, Section II.3] or
[30, Section 5.3]. As $\Psi$ satisfies the $(PS)_{c}^{G}$-condition for every
$c\in\mathbb{R}$,

$\inf_{\mathcal{U}^{G}}\Psi=\inf_{\mathcal{N}^{G}}\mathcal{J}$

is attained, i.e., the system (1.1) has a least energy $G$-invariant solution.
Moreover, since $\Psi$ is even and $\mathrm{genus}(\mathcal{U}^{G})=\infty$, a
standard variational argument shows that $\Psi$ has an unbounded sequence of
critical values. ∎

##  3 The limit profile of minimizers

We assume throughout that $1\leq\dim(Gp)<m$ for every $p\in M$. Let $\Omega$
be an open $G$-invariant subset of $M$, $H^{1}_{g,0}(\Omega)$ be the closure
of $\mathcal{C}_{c}^{\infty}(\Omega)$ in $H^{1}_{g}(M)$ and
$H^{1}_{g,0}(\Omega)^{G}$ be the space of $G$-invariant functions in
$H^{1}_{g,0}(\Omega)$. The solutions of (1.2) are the critical points of the
$\mathcal{C}^{2}$-functional $J_{\Omega}:H^{1}_{g,0}(\Omega)^{G}\to\mathbb{R}$
given by

$J_{\Omega}(u):=\frac{1}{2}\|u\|^{2}_{g}-\frac{1}{2^{*}}\int_{M}|u|^{2^{*}}\mathrm{\,d}\mu_{g}.$

The nontrivial ones belong to the Nehari manifold

$\mathcal{N}_{\Omega}^{G}:=\\{u\in H^{1}_{g,0}(\Omega)^{G}:u\neq 0\text{ and
}J_{\Omega}^{\prime}(u)u=0\\},$

which is a natural constraint for $J_{\Omega}$. So, a minimizer for
$J_{\Omega}$ over $\mathcal{N}_{\Omega}^{G}$ is a nontrivial solution of
(1.2). We call it a _least energy $G$-invariant solution_. A standard argument
using Lemma 2.1 shows that a minimizer does exist. Note that
$J_{\Omega}(u)=\frac{1}{m}\|u\|_{g}^{2}=\frac{1}{m}\int_{M}|u|^{2^{*}}\mathrm{\,d}\mu_{g}$
if $u\in\mathcal{N}_{\Omega}^{G}$, so the quantity defined in (1.3) is

$c_{\Omega}^{G}=\inf_{u\in\mathcal{N}_{\Omega}^{G}}J_{\Omega}(u).$

###### Proof of Theorem 1.2.

Let $\lambda_{n}\to-\infty$. We write $\mathcal{J}_{n}$ and
$\mathcal{N}_{n}^{G}$ for the functional and the set defined in (2.1) and
(2.2) with $\lambda_{ij}=\lambda_{n}$ for all $i\neq j$. Let

$\displaystyle\mathcal{M}_{\ell}^{G}:=$
$\displaystyle\\{(v_{1},\ldots,v_{\ell})\in\mathcal{H}^{G}:v_{i}\neq 0,\
\|v_{i}\|_{g}^{2}=\int_{M}|v_{i}|^{2^{*}}\mathrm{\,d}\mu_{g},\
v_{i}v_{j}=0\text{ on }M\text{ if }i\neq j\\},$ $\displaystyle c_{\ell}^{G}:=$
$\displaystyle\inf_{(v_{1},\ldots,v_{\ell})\in\mathcal{M}_{\ell}}\,\frac{1}{m}\sum_{i=1}^{\ell}\|v_{i}\|_{g}^{2}.$

Since $\dim(Gp)<m$ for $p\in M$, the set $\mathcal{M}_{\ell}^{G}$ is nonempty
and, so, $c_{\ell}^{G}<\infty$.

Let $\bar{u}_{n}=(u_{n,1},\ldots,u_{n,\ell})\in\mathcal{N}_{n}^{G}$ be such
that $\mathcal{J}_{n}(\bar{u}_{n})=\inf_{\mathcal{N}_{n}^{G}}\mathcal{J}_{n}$
and $u_{n,i}>0$ for every $n$ and $i$. Noting that
$\mathcal{M}_{\ell}^{G}\subset\mathcal{N}_{n}^{G}$ for every $n\in\mathbb{N}$
and recalling Lemma 2.2, we see that

$0<\inf_{\mathcal{N}_{n}^{G}}\mathcal{J}_{n}=\frac{1}{m}\sum_{i=1}^{\ell}\|u_{n,i}\|_{g}^{2}\leq
c_{\ell}^{G}<\infty\qquad\text{for every \ }n\in\mathbb{N}.$

Applying Lemma 2.1 and passing to a subsequence, we get that
$u_{n,i}\rightharpoonup u_{\infty,i}$ weakly in $H^{1}_{g}(M)^{G}$,
$u_{n,i}\to u_{\infty,i}$ strongly in $L^{2^{*}}_{g}(M)$ and $u_{n,i}\to
u_{\infty,i}$ a.e. on $M$, for each $i=1,\ldots,\ell$. Hence,
$u_{\infty,i}\geq 0$.

As $\bar{u}_{n}\in\mathcal{N}_{n}^{G}$, we have that

$\displaystyle 0$
$\displaystyle\leq\int_{M}|u_{n,j}|^{\frac{2^{*}}{2}}|u_{n,i}|^{\frac{2^{*}}{2}}\mathrm{\,d}\mu_{g}\leq\frac{\int_{M}|u_{n,i}|^{2^{*}}\mathrm{\,d}\mu_{g}}{-\lambda_{n}}\leq\frac{C}{-\lambda_{n}}\qquad\text{for
each pair \ }j\neq i,$

and, letting $n\to\infty$, we obtain

$\int_{M}|u_{\infty,j}|^{\frac{2^{*}}{2}}|u_{\infty,i}|^{\frac{2^{*}}{2}}\mathrm{\,d}\mu_{g}=0.$

So $u_{\infty,j}u_{\infty,i}=0$ a.e. on $M$ whenever $i\neq j$. We also have
that

$0<d_{0}\leq\|u_{n,i}\|_{g}^{2}\leq\int_{M}|u_{n,i}|^{2^{*}}\mathrm{\,d}\mu_{g}\qquad\text{for
every \ }n\in\mathbb{N}\text{ \ and \ }i=1,\ldots,\ell,$

with $d_{0}$ as in Lemma 2.2. Passing to the limit as $n\to\infty$, we see
that $u_{\infty,i}\neq 0$ and

$\|u_{\infty,i}\|_{g}^{2}\leq\int_{M}|u_{\infty,i}|^{2^{*}}\mathrm{\,d}\mu_{g}\qquad\text{for
every \ }i=1,\ldots,\ell.$

Hence, there exists $t_{i}\in(0,1]$ such that
$\|t_{i}u_{\infty,i}\|_{g}^{2}=\int_{M}|t_{i}u_{\infty,i}|^{2^{*}}\mathrm{\,d}\mu_{g}$.
Consequently,
$(t_{1}u_{\infty,1},\ldots,t_{\ell}u_{\infty,\ell})\in\mathcal{M}_{\ell}^{G}$
and

$c_{\ell}^{G}\leq\frac{1}{m}\sum_{i=1}^{\ell}\|t_{i}u_{\infty,i}\|_{g}^{2}\leq\frac{1}{m}\sum_{i=1}^{\ell}\|u_{\infty,i}\|_{g}^{2}\leq\frac{1}{m}\liminf_{n\to\infty}\sum_{i=1}^{\ell}\|u_{n,i}\|_{g}^{2}\leq
c_{\ell}^{G}.$

It follows that $t_{i}=1$ and $u_{n,i}\to u_{\infty,i}$ strongly in
$H^{1}_{g}(M)$. But then,
$\|u_{\infty,i}\|_{g}^{2}=\int_{M}|u_{\infty,i}|^{2^{*}}\mathrm{\,d}\mu_{g}$
and, passing to the limit in

$\displaystyle\sum_{i=1}^{\ell}\|u_{n,i}\|_{g}^{2}$
$\displaystyle=\sum_{i=1}^{\ell}|u_{n,i}|_{g,2^{*}}^{2^{*}}+\mathop{\sum_{i,j=1}^{\ell}}_{j\neq
i}\int_{M}\lambda_{n}|u_{n,j}|^{\frac{2^{*}}{2}}|u_{n,i}|^{\frac{2^{*}}{2}},$

we obtain

$\lim_{n\to\infty}\int_{M}\lambda_{n}|u_{\infty,j}|^{\frac{2^{*}}{2}}|u_{\infty,i}|^{\frac{2^{*}}{2}}=0\qquad\text{for
every \ }i\neq j.$

Moreover, $(u_{\infty,1},\ldots,u_{\infty,\ell})\in\mathcal{M}_{\ell}^{G}$ and

$\frac{1}{m}\sum_{i=1}^{\ell}\|u_{\infty,i}\|_{g}^{2}=c_{\ell}^{G}.$

We have now all assumptions needed to apply [9, Lemmas 4.3 and 4.4] and
conclude that $(u_{n,i})$ is uniformly bounded in the Hölder norm, i.e., for
any $\alpha\in(0,1)$ there exists $C_{\alpha}>0$ such that

$\|u_{n,i}\|_{\mathcal{C}^{0,\alpha}(M)}\leq C_{\alpha}\qquad\text{for every \
}n\in\mathbb{N}\text{ \ and \ }i=1,\ldots,\ell.$

As a consequence, $u_{n,i}\to u_{\infty,i}$ in $\mathcal{C}^{0,\alpha}(M)$ for
every $i$. In particular, $u_{\infty,i}$ is continuous and $G$-invariant, so
the set $\Omega_{i}:=\\{p\in M:u_{\infty,i}(p)>0\\}$ is open and
$G$-invariant. Since $u_{\infty,i}u_{\infty,j}=0$ if $i\neq j$, we have that
$\Omega_{i}\cap\Omega_{j}=\emptyset$. This shows that
$\\{\Omega_{1},\ldots,\Omega_{\ell}\\}\in\mathcal{P}_{\ell}^{G}$.

We claim that $u_{\infty,i}$ is a least energy $G$-invariant solution to (1.2)
in $\Omega_{i}$ for all $i$. Otherwise,
$J_{\Omega_{i}}(u_{\infty,i})>c_{\Omega_{i}}^{G}$ for some $i$ and there would
exist $v_{i}\in\mathcal{N}_{\Omega_{i}}^{G}$ with
$c_{\Omega_{i}}^{G}<J_{\Omega_{i}}(v_{i})<J_{\Omega_{i}}(u_{\infty,i})$. But
then, setting $v_{j}:=u_{\infty,j}$ for $j\neq i$, we would have that
$(v_{1},\ldots,v_{\ell})\in\mathcal{M}_{\ell}^{G}$ and

$\frac{1}{m}\sum_{i=1}^{\ell}\|v_{i}\|_{g}^{2}<\frac{1}{m}\sum_{i=1}^{\ell}\|u_{\infty,i}\|_{g}^{2}=c_{\ell}^{G},$

contradicting the definition of $c_{\ell}^{G}$. Hence,
$J_{\Omega_{i}}(u_{\infty,i})=c_{\Omega_{i}}^{G}$ for all $i=1,\ldots,\ell$,
as claimed. A similar argument shows that $\Omega_{i}$ is connected.

As a consequence, if
$\\{\Theta_{1},\ldots,\Theta_{\ell}\\}\in\mathcal{P}_{\ell}^{G}$, taking
$w_{i}\in\mathcal{N}_{\Theta_{i}}^{G}$ with
$J_{\Theta_{i}}(w_{i})=c_{\Theta_{i}}^{G}$, we have that
$(w_{1},\ldots,w_{\ell})\in\mathcal{M}_{\ell}^{G}$ and, therefore,

$\sum_{i=1}^{\ell}c_{\Omega_{i}}^{G}=\frac{1}{m}\sum_{i=1}^{\ell}\|u_{\infty,i}\|_{g}^{2}=c_{\ell}^{G}\leq\frac{1}{m}\sum_{i=1}^{\ell}\|w_{i}\|_{g}^{2}=\sum_{i=1}^{\ell}c_{\Theta_{i}}^{G}.$

This shows that $\\{\Omega_{1},\ldots,\Omega_{\ell}\\}$ is an optimal
$(G,\ell)$-partition, and completes the proof of statements $(i)$ and $(iii)$.

Statements $(ii)$ and $(iv)$ are local. In local coordinates the system (1.1)
with $\lambda_{ij}$ replaced by $\lambda_{n}$ becomes

$-\mathrm{div}(A(x)\nabla
v_{i})=f_{i}(x,v_{i})+a(x)\sum_{\begin{subarray}{c}j=1\\\ j\neq
i\end{subarray}}^{\ell}\lambda_{n}|v_{j}|^{\frac{2^{*}}{2}}|v_{i}|^{\frac{2^{*}}{2}-2}v_{i},\qquad
x\in\Omega,$

where $\Omega\subset\mathbb{R}^{m}$ is open and bounded,
$a(x):=\sqrt{|g(x)|}$, $A(x):=\sqrt{|g(x)|}(g^{kl}(x))$,
$f_{i}(x,s):=a(x)(|s|^{2^{*}-2}s-\kappa_{m}S_{g}(x)s)$ and, as usual, $|g|$
denotes the determinant of the metric $g=(g_{kl})$ in local coordinates and
$(g^{kl})$ its inverse. Applying [9, Theorem C.1] we conclude that $(ii)$ and
$(iv)$ are true locally on $M$, hence also globally.

A standard argument yields the proof of statement $(v)$. Namely, the
$G$-invariant sign-changing solutions to the Yamabe equation (1.4) belong to
the set

$\mathcal{E}_{M}^{G}:=\\{u\in\mathcal{N}_{M}^{G}:u^{+}\in\mathcal{N}_{M}^{G}\text{
\ and \ }u^{-}\in\mathcal{N}_{M}^{G}\\},$

where $u^{+}:=\max\\{u,0\\}\neq 0$ and $u^{-}:=\min\\{u,0\\}\neq 0$. Moreover,
as shown in [4, Lemma 2.6], any minimizer of $J_{M}$ on $\mathcal{E}_{M}^{G}$
is a $G$-invariant sign-changing solution of (1.4). For every
$u\in\mathcal{E}_{M}^{G}$, we have that $(u^{+},u^{-})\in\mathcal{M}_{2}^{G}$
and $J_{M}(u)=\frac{1}{m}(\|u^{+}\|_{g}^{2}+\|u^{-}\|_{g}^{2})$. Therefore,

$\inf_{\mathcal{E}_{M}^{G}}J_{M}\geq
c_{2}^{G}=\frac{1}{m}(\|u_{\infty,1}\|_{g}^{2}+\|u_{\infty,2}\|_{g}^{2}).$

As $u_{\infty,1}-u_{\infty,2}\in\mathcal{E}_{M}^{G}$, it is a minimizer of
$J_{M}$ on $\mathcal{E}_{M}^{G}$. This completes the proof. ∎

##  4 The Yamabe system on the standard sphere

In this section we give an account of some known results for the system (1.1)
on the standard sphere.

The stereographic projection
$\sigma:\mathbb{S}^{m}\smallsetminus\\{N\\}\to\mathbb{R}^{m}$ from the north
pole $N$ is a conformal diffeomorphism. The conformal invariance of the
operator $\mathscr{L}_{\bar{g}}$ (see [23, Proposition 6.1.1]) allows to
establish a one-to-one correspondence between solutions to the system (1.1) on
the standard sphere $\mathbb{S}^{m}$ and solutions to the system

$\begin{cases}-\Delta
v_{i}=|v_{i}|^{2^{*}-2}v_{i}+\sum\limits_{\begin{subarray}{c}j=1\\\ j\neq
i\end{subarray}}^{\ell}\lambda_{ij}|v_{j}|^{\frac{2^{*}}{2}}|v_{i}|^{\frac{2^{*}}{2}}v_{i}\\\
v_{i}\in D^{1,2}(\mathbb{R}^{m}),\quad v_{i}\neq 0,\qquad
i=1,\ldots,\ell,\end{cases}$ (4.1)

where $D^{1,2}(\mathbb{R}^{m}):=\\{v\in L^{2^{*}}(\mathbb{R}^{m}):\nabla v\in
L^{2}(\mathbb{R}^{m},\mathbb{R}^{m})\\}$; see [7, 10].

Let us first consider the case of the single equation

$-\Delta u=|u|^{2^{*}-2}u,\qquad u\in D^{1,2}(\mathbb{R}^{m}).$ (4.2)

It is well known [3, 25, 29] that all positive solutions to (4.2) are the so-
called standard bubbles

$U_{\delta,y}(x):=\mathfrak{c}_{m}\frac{\delta^{\frac{m-2}{2}}}{(\delta^{2}+|x-y|^{2})^{\frac{m-2}{2}}},\qquad
x,y\in\mathbb{R}^{m},\ \delta>0.$

The existence of sign-changing solutions was first established by W.Y. Ding in
[17]. He considered solutions that are invariant under the action of the group
$\Gamma:=O(n_{1})\times O(n_{2})$ with $n_{1},n_{2}\geq 2$ and
$n_{1}+n_{2}=m+1$, acting on
$\mathbb{S}^{m}\subset\mathbb{R}^{m+1}\equiv\mathbb{R}^{n_{1}}\times\mathbb{R}^{n_{2}}$
in the obvious way. Using a variational approach he proved the existence of
infinitely many $\Gamma$-invariant sign-changing solutions to the Yamabe
equation (1.4) on $\mathbb{S}^{m}$ or, equivalently, to the equation (4.2) on
$\mathbb{R}^{m}$.

In [15, 16], del Pino, Musso, Pacard and Pistoia exploited the symmetries of
the sphere to build solutions to (4.2) which have large energy and concentrate
along some special submanifolds of $\mathbb{S}^{m}$. In particular, for $m\geq
4$ they obtained sequences of solutions to the Yamabe equation whose energy
concentrates along one great circle or finitely many great circles that are
linked to each other (corresponding to Hopf links embedded in
$\mathbb{S}^{3}\times\\{0\\}\subset\mathbb{S}^{m}$), and for $m\geq 5$ they
also obtained sequences of solutions whose energy concentrates along a two-
dimensional Clifford torus in
$\mathbb{S}^{3}\times\\{0\\}\subset\mathbb{S}^{m}$. These solutions are built
via gluing techniques (e.g., a Ljapunov-Schmidt procedure) and can be
described as the superposition of the constant solution to (1.4) with a large
number of copies of negative solutions of (1.4) which are highly concentrated
at points evenly arranged along some special submanifolds of the sphere.

Recently, Fernández and Petean in [19] established the existence of solutions
to the Yamabe problem (1.4) on $\mathbb{S}^{m}$ with precisely $\ell$ nodal
domains which are isoparametric hypersurfaces, for every $\ell\geq 2$, using
ODE techniques.

Regarding the competitive system (where all $\lambda_{ij}<0$), Guo, Li and Wei
[21] established the existence of solutions for (4.1) with $\ell=2$ and
$\lambda_{12}<0$ in $\mathbb{R}^{3}$, and using the approach developed in [15]
they built a sequence of positive nonradial solutions whose first component
looks like the constant function and the second component resembles a large
number of copies of positive solutions of (1.4) concentrated at points that
are placed along a circle. The argument of their proof, which relies on the
Ljapunov-Schmidt procedure, cannot be extended to higher dimensions because
the coupling terms have linear (if $m=4$) or sublinear (if $m\geq 5$) growth.

Following the variational approach presented in the previous sections,
successively Clapp and Pistoia [7], Clapp and Szulkin [11] and Clapp, Saldaña
and Szulkin [10] found $\Gamma$-invariant solutions to the competitive Yamabe
system (1.1) on $\mathbb{S}^{m}$ for the groups considered by Ding, and
described the limit profile of least energy solutions as
$\lambda_{ij}\to\\-\infty$. In this particular case, Theorems 1.1 and 1.2 were
proved in [7, 11] and [7, 10] respectively. In addition, a more accurate
description of the optimal $(\Gamma,\ell)$-partition of $\mathbb{S}^{m}$ is
provided in [7, 10]. Note that the $\Gamma$-orbit of a point
$p\in\mathbb{S}^{m}$ is diffeomorphic to either $\mathbb{S}^{n_{1}-1}$, or
$\mathbb{S}^{n_{1}-1}\times\mathbb{S}^{n_{2}-1}$, or $\mathbb{S}^{n_{2}-1}$,
and that the map $q:\mathbb{S}^{m}\to[0,\pi]$ given by

$q(x,y):=\arccos(|x|^{2}-|y|^{2}),\quad\text{where \ }x\in\mathbb{R}^{n_{1}},\
y\in\mathbb{R}^{n_{2}},$ (4.3)

is a quotient map identifying each $\Gamma$-orbit in $\mathbb{S}^{m}$ to a
single point. So $q$ maps a $(\Gamma,\ell)$-partition of $\mathbb{S}^{m}$ onto
a partition of $[0,\pi]$ by relatively open subintervals. This last partition
can be ordered. Taking advantage of this fact, the following result was proved
in [7, 10].

######  Theorem 4.1.

Let $\Gamma:=O(n_{1})\times O(n_{2})$ with $n_{1},n_{2}\geq 2$ and
$n_{1}+n_{2}=m+1$ and let $M:=\mathbb{S}^{m}$. Then, the optimal
$(\Gamma,\ell)$-partition
$\\{\Omega_{1},\ldots,\Omega_{\ell}\\}\in\mathcal{P}_{\ell}^{\Gamma}$ of
$\mathbb{S}^{m}$ given by Theorem 1.2 has the following properties:
$\Omega_{1},\ldots,\Omega_{\ell}$ are smooth and connected,
$\overline{\Omega_{1}\cup\cdots\cup\Omega_{\ell}}=\mathbb{S}^{m}$ and, after
reordering,

  * •

$\Omega_{1}\cong\mathbb{S}^{n_{1}-1}\times\mathbb{B}^{n_{2}}$,
$\Omega_{i}\cong\mathbb{S}^{n_{1}-1}\times\mathbb{S}^{n_{2}-1}\times(0,1)$ if
$i=2,\ldots,\ell-1$, and
$\Omega_{\ell}\cong\mathbb{B}^{n}_{1}\times\mathbb{S}^{n_{2}-1}$,

  * •

$\overline{\Omega}_{i}\cap\overline{\Omega}_{i+1}\cong\mathbb{S}^{n_{1}-1}\times\mathbb{S}^{n_{2}-1}$
and $\overline{\Omega}_{i}\cap\overline{\Omega}_{j}=\emptyset$ if $|j-i|\geq
2$,

  * •

the function

$u:=\sum_{i=1}^{\ell}(-1)^{i-1}u_{\infty,i}$

is a $\Gamma$-invariant sign-changing solution to the Yamabe problem (1.4) on
$\mathbb{S}^{m}$ with precisely $\ell$ nodal domains, and $u$ has least energy
among all such solutions.

As shown by Theorems 1.2 and 4.1 the competitive Yamabe system (1.1) gives
rise to one or even multiple solutions to the Yamabe equation (1.4). The
opposite question has been also considered in the literature. Let us call a
solution $\mathbf{u}=(u_{1},\ldots,u_{\ell})$ to the Yamabe system (4.1)
_fully synchronized_ if there exist $c_{i}\neq 0$ and a nontrivial solution
$u$ to the single equation (4.2) such that $u_{i}=c_{i}u$ for all
$i=1,\ldots,\ell$. It is readily seen that $(c_{1}u,\ldots,c_{\ell}u)$ solves
(4.1) for any solution $u$ of (4.2) iff
$\mathbf{c}=(c_{1},\ldots,c_{\ell})\in\mathbb{R}^{\ell}$ solves the algebraic
system

$c_{i}=|c_{i}|^{2^{*}-2}c_{i}+\sum_{\begin{subarray}{c}j=1\\\
j\not=1\end{subarray}}^{\ell}\lambda_{ij}|c_{j}|^{\frac{2^{*}}{2}}|c_{i}|^{\frac{2^{*}}{2}-2}c_{i},\quad
c_{i}>0,\quad\text{for every \ }i=1,\dots,\ell.$ (4.4)

There are several results concerning the solvability of (4.4). The easiest
case is when $m=4$ (i.e., $2^{*}=4$) and $\ell=2$. Indeed, a straightforward
computation shows that a solution to (4.4) exists if and only if
$\lambda_{12}>-1$ and $\lambda_{12}\not=1$. Bartsch proved in [2, Proposition
2.1] that, if $2^{*}=4$ and $\ell\geq 2$, a fully synchronized solution to
(4.4) exists when $\lambda_{ij}:=\lambda\not=1$ for all $i\not=j$ and
$\lambda>\overline{\lambda}$ for some $\overline{\lambda}<0$ , while Chen and
Zou [14, Theorem 1.1] showed that if $m\geq 5$ (i.e., $2^{*}<4$) and $\ell=2$
a fully synchronized solution to (4.4) always exists provided
$\lambda_{12}>0$. Recently, Clapp and Pistoia complemented these results in
[8] showing that, if the system is purely cooperative (i.e., $\lambda_{ij}\geq
0$ for all $i,j=1,\ldots,\ell$, $i\neq j$), there exists a solution to (4.4).
We recall that in the purely cooperative case every positive solution of (4.1)
with $\ell=2$ is fully syncronized, as shown by Guo and Liu in [22]. On the
other hand, it is shown in [7, 11] that there exists $\lambda^{*}<0$ such that
the system (4.1) does not have a fully synchronized solution if
$\lambda_{ij}<\lambda^{*}$ for all pairs $i\neq j$.

System (4.1) with mixed couplings (i.e., $\lambda_{ij}$ can be positive or
negative) has been recently studied by Clapp and Pistoia [8]. Let
$\lambda_{ii}=1$ and assume the matrix $(\lambda_{ij})$ is symmetric and
admits a block decomposition as follows: For some $1<q<\ell$ there exist
$0=\ell_{0}<\ell_{1}<\dots<\ell_{q-1}<\ell_{q}=\ell$ such that, if we set

$\displaystyle I_{h}:=\\{i\in\\{1,\dots,\ell\\}:\ell_{h-1}<i\leq\ell_{h}\\},$
$\displaystyle\mathcal{I}_{h}:=I_{h}\times
I_{h},\qquad\mathcal{K}_{h}:=\big{\\{}(i,j)\in I_{h}\times
I_{k}:k\in\\{1,\ldots,q\\}\smallsetminus\\{h\\}\big{\\}},$

then

$\lambda_{ij}>0\ \text{ if }\ (i,j)\in\mathcal{I}_{h}\quad\text{ and
}\quad\lambda_{ij}\leq 0\ \text{if}\ (i,j)\in\mathcal{K}_{h},\quad
h=1,\ldots,q.$

According to the above decomposition, we shall write a solution
$\mathbf{u}=(u_{1},\ldots,u_{\ell})$ to (4.1) in block-form as

$\mathbf{u}=(\bar{u}_{1},\ldots,\bar{u}_{q})\qquad\text{with \
}\bar{u}_{h}=(u_{\ell_{h-1}+1},\ldots,u_{\ell_{h}}).$

$\mathbf{u}$ is called _fully nontrivial_ if every component $u_{i}$ is
different from zero. In [8] it is proved that if, either $m\geq 5$, or $m=4$
and $\lambda_{ij}=:b_{h}>1$ for all $i,j\in I_{h}$ with $i\neq j$, the system
(4.1) has a fully nontrivial solution if
$\max_{(i,j)\in\mathcal{K}_{h}}|\lambda_{ij}|<\Lambda$ for some $\Lambda>0$.
This solution is invariant under the conformal action on $\mathbb{R}^{m}$ of
the group $\Gamma$ defined above.

Finally, we would like to mention a couple of results, one of them by Grossi,
Gladiali and Troestler [20] where they give sufficient conditions on the
matrix $(\lambda_{ij})$ to ensure the existence of solutions bifurcating from
the bubble of the critical Sobolev equation, and another one by Druet and
Hebey [18] where they study the stability of solutions to (4.1) under linear
perturbation.

## References

  * [1] Ammann, Bernd; Humbert, Emmanuel: The second Yamabe invariant. J. Funct. Anal. 235 (2006), no. 2, 377–412.
  * [2] Bartsch, Thomas: Bifurcation in a multicomponent system of nonlinear Schrödinger equations. J. Fixed Point Theory Appl. 13 (2013), no. 1, 37–50.
  * [3] Aubin, Thierry: Problèmes isopérimétriques et espaces de Sobolev, J. Differ. Geom. 11 (4) (1976) 573–598.
  * [4] Castro, Alfonso; Cossio, Jorge; Neuberger, John M.: A sign-changing solution for a superlinear Dirichlet problem. Rocky Mountain J. Math. 27 (1997), no. 4, 1041–1053.
  * [5] Chang, Shu-Ming; Lin, Chang-Shou; Lin, Tai-Chia, Lin, Wen-Wei: Segregated nodal domains of two-dimensional multispecies Bose-Einstein condensates. Phys. D 196 (2004), 341–361.
  * [6] Clapp, Mónica; Fernández, Juan Carlos: Multiplicity of nodal solutions to the Yamabe problem. Calc. Var. Partial Differential Equations 56 (2017), no. 5, Paper No. 145, 22 pp.
  * [7] Clapp, Mónica; Pistoia, Angela: Existence and phase separation of entire solutions to a pure critical competitive elliptic system. Calc. Var. Partial Differential Equations 57 (2018), no. 1, Paper No. 23, 20 pp.
  * [8] Clapp, Mónica; Pistoia, Angela: Fully nontrivial solutions to elliptic systems with mixed couplings (2021) ArXiv:2106.01637
  * [9] Clapp, Mónica; Pistoia, Angela; Tavares, Hugo: Yamabe systems, optimal partitions and nodal solutions to the Yamabe equation. Preprint 2021, arXiv:2106.00579.
  * [10] Clapp, Mónica; Saldaña, Alberto; Szulkin, Andrzej: Phase separation, optimal partitions and nodal solutions to the Yamabe equation on the sphere. International Mathematics Research Notices 2021 (2021), no. 5, 3633–3652.
  * [11] Clapp, Mónica; Szulkin, Andrzej: A simple variational approach to weakly coupled competitive elliptic systems. Nonlinear Differential Equations and Applications NoDEA 26:26 (2019), 21 pp.
  * [12] Conti, Monica; Terracini, Susanna; Verzini, Gianmaria: Nehari’s problem and competing species systems. Ann. Inst. H. Poincaré Anal. Non Linéaire 19 (2002), no. 6, 871–888.
  * [13] Conti, Monica; Terracini, Susanna; Verzini, Gianmaria: A variational problem for the spatial segregation of reaction-diffusion systems. Indiana Univ. Math. J. 54 (2005), no. 3, 779–815.
  * [14] Chen, Zhijie; Zou, Wenming: Positive least energy solutions and phase separation for coupled Schrödinger equations with critical exponent: higher dimensional case. Calc. Var. Partial Differential Equations 52 (2015), no. 1-2, 423-467.
  * [15] del Pino, Manuel; Musso, Monica; Pacard, Frank; Pistoia, Angela: Large energy entire solutions for the Yamabe equation. J. Differential Equations 251 (9) (2011) 2568–2597.
  * [16] del Pino, Manuel; Musso, Monica; Pacard, Frank; Pistoia, Angela: Torus action on Sn and sign-changing solutions for conformally invariant equations. Ann. Sc. Norm. Super. Pisa Cl. Sci. (5) 12 (1) (2013) 209–237.
  * [17] Ding, Wei Yue: On a conformally invariant elliptic equation on $R^{n}$. Comm. Math. Phys. 107 (1986), no. 2, 331–335.
  * [18] Druet, Olivier; Hebey, Emmanuel: Stability for strongly coupled critical elliptic systems in a fully inhomogeneous medium. Anal. PDE 2 (2009), no. 3, 305–359.
  * [19] Fernández, Juan Carlos; Petean, Jimmy: Low energy nodal solutions to the Yamabe equation. J. Differential Equations 268 (2020), no. 11, 6576–6597.
  * [20] Gladiali, Francesca; Grossi, Massimo; Troestler, Christophe: A non-variational system involving the critical Sobolev exponent. The radial case. J. Anal. Math. 138 (2019), no. 2, 643–671.
  * [21] Guo, Yuxia; Li, Bo; Wei, Juncheng: Entire nonradial solutions for non-cooperative coupled elliptic system with critical exponents in $\mathbb{R}^{3}$. J. Differential Equations 256 (2014), no. 10, 3463–3495.
  * [22] Guo, Yuxia; Liu, Jiaquan: Liouville type theorems for positive solutions of elliptic system in $\mathbb{R}^{n}$. Comm. Partial Differential Equations 33 (2008), no. 1–3, 263–284.
  * [23] Hebey, E.: Introduction à l’analyse non linéaire sur les variétés. Diderot, Paris, 1997.
  * [24] Hebey, Emmanuel; Vaugon, Michel: Sobolev spaces in the presence of symmetries. J. Math. Pures Appl. (9) 76 (1997), no. 10, 859–881.
  * [25] Obata, Morio: The conjectures on conformal transformations of Riemannian manifolds. J. Differential Geometry 6 (1971/72), 247–258.
  * [26] Palais, Richard S.: The principle of symmetric criticality. Comm. Math. Phys. 69 (1979), no. 1, 19–30.
  * [27] Struwe, Michael: Variational methods. Applications to nonlinear partial differential equations and Hamiltonian systems. Second edition. Ergebnisse der Mathematik und ihrer Grenzgebiete, 34. Springer-Verlag, Berlin, 1996.
  * [28] Soave, Nicola; Tavares, Hugo; Terracini, Susanna; Zilio, Alessandro: Hölder bounds and regularity of emerging free boundaries for strongly competing Schrödinger equations with nontrivial grouping. Nonlinear Anal. 138 (2016), 388–427.
  * [29] Talenti, Giorgio: Best constant in Sobolev inequality. Ann. Mat. Pura Appl. (4) 110 (1976), 353–372.
  * [30] Willem, Michel: Minimax theorems. Progress in Nonlinear Differential Equations and their Applications, 24. Birkhäuser Boston, Inc., Boston, MA, 1996.

Mónica Clapp
Instituto de Matemáticas

Universidad Nacional Autónoma de México

Circuito Exterior, Ciudad Universitaria

04510 Coyoacán, Ciudad de México, Mexico

<EMAIL_ADDRESS>
Angela Pistoia
Dipartimento di Metodi e Modelli Matematici

La Sapienza Università di Roma

Via Antonio Scarpa 16

00161 Roma, Italy

<EMAIL_ADDRESS>


# Reducing computational effort in topology optimization considering the
deformation in additive manufacturing.

Takao Miki<EMAIL_ADDRESS>Osaka Research Institute of Industrial Science and
Technology, 7-1, Ayumino-2, Izumi-city, Osaka, 594-1157, Japan

###### Abstract

Integrating topology optimization and additive manufacturing (AM) technology
can facilitate innovative product development. However, laser powder bed
fusion, which is the predominant method in metal AM, can lead to issues such
as residual stress and deformation. Recently, topology optimization methods
considering these stresses and deformations have been proposed; however, they
suffer from challenges caused by an increased computational cost. In this
study, we propose a method for reducing computational cost in topology
optimization considering the deformation in AM. An inherent strain method-
based analytical model is presented for simulating the residual stress and
deformation in the AM process. Subsequently, a constraint condition to
suppress the deformation is formulated, and a method to reduce the
computational cost of the adjoint analysis in deriving sensitivity is
proposed. The minimum mean compliance problem considering AM deformation and
self-support constraints can then be incorporated into the level set-based
topology optimization framework. Finally, numerical examples are presented for
validating the effectiveness of the proposed topology optimization method.

###### keywords:

Topology optimization, Level set method, Laser powder bed fusion (metal
additive manufacturing), Design for additive manufacturing, Computational cost
, Inherent strain method

††journal: Journal of LaTeX Templates

##  1 Introduction

Additive manufacturing (AM) processes can fabricate complex geometric parts by
repeating the process of building up materials and be expected to further
improve product performance. In recent years, several research studies have
focused on innovative manufacturing combined with topology optimization, which
is a powerful design tool that can exploit the advantages of the features of
this manufacturing technique.

A prevalent method in metal AM is laser powder bed fusion (LPBF), which can
construct three-dimensional objects by sequentially melting and solidifying
metal powder using a laser in a layer-by-layer manner. Despite its widespread
use, LPBF suffers from several manufacturing challenges, including overhang
limitation, residual stress, distortion, and overheating. Support structures
are often required for mitigating issues compromising dimensional accuracy and
inducing manufacturing defects. However, these structures are not ideal
because they elevate manufacturing and post-processing costs because of their
removal.

Topology optimization methods that can address these manufacturing challenges
of AM have been developed actively. Topology optimization methods that can
incorporate overhang constraints have received considerable attention within
these challenges [1, 2, 3, 4, 5, 6, 7, 8, 9]. A common strategy involves
integrating geometrical constraints for creating self-supporting structures,
minimizing their need for support structures and associated costs. However,
this strategy cannot adequately address physical issues such as residual
stress and deformation.

Several topology optimization approaches have been explored to consider
residual stress, deformation, and overheating [10, 11, 12, 13, 14, 15, 16,
17]. These methods can circumvent physical challenges by integrating AM
process models, and they can be categorized into the thermo-mechanical
analysis [18, 19, 20, 21, 22, 23, 24] and inherent strain methods [25, 26, 27,
28, 29]. The former employs a coupled transient thermal and mechanical
analysis for forecasting component overheating, residual stress, and
deformation. The latter predicts residual stress and deformation using an
inherent strain derived from the thermo-mechanical analysis or experimental
data; however, it does not consider overheating caused by the absence of
thermal information. The main advantage lies in its lower computational
demand.

Some topology optimization methods consider multiple manufacturability
factors, thereby merging geometrical constraints with AM process models. Miki
[30] proposed a level-set-based topology optimization method that incorporates
overhang and deformation constraints. Xu et al. [31] proposed a density-based
topology optimization method with overhang and residual stress constraints.
Both approaches utilize a simplified inherent strain method that excludes
plastic deformation considerations; however, they are computationally
intensive because of iterative calculations required to optimize residual
stress or deformation in the AM process, thereby necessitating layer-by-layer
mechanical and adjoint analyses.

A potential solution for reducing computational burden involves increasing the
macro layer size in the AM process model, decreasing the number of layers.
Although this approach decreases the computational load of the mechanical and
adjoint analyses, it also affects the accuracy of residual stress and
deformation predictions. To overcome this limitation, we propose a method for
reducing computational costs while preserving prediction accuracy. This is
achieved by streamlining the adjoint analysis to a single calculation,
focusing on the adjoint problem associated with the final state of the
construction process of AM, instead of each intermediate state.

The remainder of this paper is organized as follows. Section 2 describes the
level-set-based topology optimization method. In Section 3, we present the AM
process model based on the inherent strain method and a physical model
representing the final state of the AM building process. Subsequently, the
formulation of the constraint function for physical challenges in AM is
introduced. In Section 4, the self-support constraint based on the combination
of geometrical and thermal constraints are formulated. Then, we formulate an
optimization problem that considers multiple challenges of the AM process in
Section 5, and construct an optimization algorithm using the finite element
method (FEM) in Section 6. Section 7 presents the 2D and 3D numerical examples
used for validating the utility of the proposed method. Finally, Section 8
provides the conclusions.

##  2 Level-set-based topology optimization

The basic concept of topology optimization involves replacing the structural
optimization problem with a material distribution problem, wherein the
presence or absence of a material is determined within a fixed design domain
$D\subset\mathbb{R}^{N}(withN=2\,or\,3)$. The fixed design domain $D$ is
differentiated into the material domain $\mathit{\Omega}$ and void domain
$D\backslash\mathit{\Omega}$ using a characteristic function. The level set
method expresses the structural boundary $\partial\mathit{\Omega}$ between the
material and void domains using an iso-surface of a scalar function called the
level set function. In this study, topology optimization [32] performed using
the reaction-diffusion equation is employed to update the level set function,
which is defined as follows:

$\left\\{\begin{array}[]{ l l }{0<\phi(\bm{x})\leq 1}&{\text{ for
}\bm{x}\in\mathit{\Omega}}\\\ {\phi(\bm{x})=0}&{\text{ for
}\bm{x}\in\partial\mathit{\Omega}}\\\ {-1\leq\phi(\bm{x})<0}&{\text{ for
}\bm{x}\in D\backslash\mathit{\Omega}}.\end{array}\right.$ (1)

The characteristic function is represented using the level-set function $\phi$
as follows:

$\chi_{\phi}=\left\\{\begin{array}[]{ll}1&\text{ for }\phi(x)\geq 0\\\
0&\text{ for }\phi(x)<0,\end{array}\right.$ (2)

based on Eq. 2, the level-set-based optimization problem is formulated as
follows:

$\begin{split}&\inf_{\phi}\quad J\left(u,\chi_{\phi}\right)\\\ &\text{ subject
to: }\quad\text{ governing equations for}\quad u\end{split}$ (3)

where $J$ and $u$ represent an objective function and a state variable,
respectively. Then, the optimal material distribution $\chi_{\phi}$, i.e., the
distribution of the level-set function $\phi$, is determined by solving a
time-evolution equation, where a fictitious time $s$ is introduced as follows:

$\displaystyle\frac{\partial\phi(s)}{\partial
s}=-K(\tilde{C}J^{\prime}-\tau\nabla^{2}\phi),$ (4)
$\displaystyle\tilde{C}:=\dfrac{C\int_{D}\mathrm{~{}d}\mathit{\Omega}}{\int_{D}|J^{\prime}|\mathrm{~{}d}\mathit{\Omega}}.$
(5)

where $K\in\mathbb{R}^{+}$, $C\in\mathbb{R}^{+}$, $J^{\prime}$, and
$\tau\in\mathbb{R}^{+}$ represent a proportional coefficient, normalization
coefficient, design sensitivity of the target optimization problem, and a
regularization parameter, respectively. The design sensitivity is derived
based on an adjoint variable method. In this study, we set $K=1.0$, $C=0.8$,
and $\tau=5.0\times 10^{-4}$.

##  3 AM process model and deformation constraint

The inherent strain method is widely used to predict residual stress and
deformation in welding. The AM process can be considered an iterative welding
process, and therefore, the inherent strain method is applicable. The residual
stress is attributed to inelastic strains such as thermal, plastic, and phase
transformation strains, which are generated during the local melting and
solidification processes of the metal powder. The sum of these inelastic
strains is called the inherent strain $\bm{\varepsilon}^{inh}$, and the total
strain
$\bm{\varepsilon}(\bm{u})\coloneqq\frac{1}{2}(\nabla\bm{u}+(\nabla\bm{u})^{\top})$
is composed of elastic strain $\bm{\varepsilon}^{el}$ and inherent strain
$\bm{\varepsilon}^{inh}$ as follows:

$\bm{\varepsilon}(\bm{u})=\bm{\varepsilon}^{el}+\bm{\varepsilon}^{inh}.$ (6)

We consider a part $\mathit{\Omega}$ divided into $m$ macro-layers with a
fixed layer thickness in the building direction to apply the inherent strain
method to the AM building process shown in Fig. 1.

Figure 1: Domains and boundaries of the mechanical model in the building
process.

$\mathit{\Omega_{p}}$ represents the printed domain and
$\mathit{\Omega_{inh}}$ represents the layer $\mathit{\Omega_{i}}$ at the top
of the printed domain $\mathit{\Omega_{p}}$, where the inherent strain is
imposed. The printed domain $\mathit{\Omega_{p}}$ is defined by domain
$\mathit{\Omega}_{i}$, $1\leq i\leq m$ as follows:

$\mathit{\Omega_{p}}=\mathit{\Omega_{1}}\cup{\ldots}\cup{\mathit{\Omega_{i}}}.$
(7)

This study employs a linear inherent strain problem that does not consider
plastic deformation [12, 31, 33]. $\mathit{\Omega_{p}}$ is occupied by an
isotropic elastic material. The bottom surface $\mathit{\Gamma_{u}}$ of the
part is clamped and $\mathit{\Omega_{inh}}$ is loaded with the inherent
strain. The displacement field $\bm{u}_{i}\in H^{1}(\mathit{\Omega}_{p})^{N}$
is represented by the following governing equation:

$\displaystyle\begin{cases}\hskip 5.69054pt-\mbox{\rm
div}(\bm{\sigma}_{i})=0\hskip 14.22636pt&\text{in}\hskip
5.69054pt\mathit{\Omega_{p}}\\\ \hskip
5.69054pt\bm{\sigma}_{i}=\mathbb{D}(\bm{\varepsilon}(\bm{u}_{i})-\bm{\varepsilon}^{inh}),\hskip
14.22636pt&\\\ \hskip 5.69054pt\bm{u}_{i}=\bm{0}\hskip
14.22636pt&\text{on}\hskip 5.69054pt\mathit{\Gamma_{u}}\\\ \hskip
5.69054pt-\bm{\sigma}_{i}\cdot\bm{n}=\bm{0}\hskip 14.22636pt&\text{on}\hskip
5.69054pt\partial\mathit{\Omega_{p}}\setminus\mathit{\Gamma_{u}},\end{cases}$
(8)

for all indices $i=1,2,\ldots,m$, where $\bm{\sigma}$ and $\mathbb{D}$
represent the stress tensor and isotropic elasticity tensor, respectively. The
inherent strain $\bm{\varepsilon}^{inh}$ in $\mathit{\Omega_{inh}}$ is defined
as follows:

$\bm{\varepsilon}^{inh}(\bm{x})=\left\\{\begin{array}[]{ l l
}{\bm{\varepsilon}^{inh}}&{\text{ for }\bm{x}\in\mathit{\Omega_{inh}}},\\\
{\bm{0}}&{\text{ otherwise }},\end{array}\right.$ (9)

The inherent strain vector $\bm{\varepsilon}^{inh}$ has six components
$[\varepsilon^{inh}_{x}\quad\varepsilon^{inh}_{y}\quad\varepsilon^{inh}_{z}\quad\varepsilon^{inh}_{xy}\quad\varepsilon^{inh}_{yz}\quad\varepsilon^{inh}_{zx}]^{\mathsf{T}}$.
In this study, the in-plane strain component
$\varepsilon^{inh}_{x},\varepsilon^{inh}_{y}$ is assumed to be constant, and
the building direction $\varepsilon^{inh}_{z}$ and shear strain components
$\varepsilon^{inh}_{xy},\varepsilon^{inh}_{yz},\varepsilon^{inh}_{zx}$ are
omitted [34, 26]. Part-scale residual stress and distortion can be calculated
by applying inherent strain $\bm{\varepsilon}^{inh}$ to each domain
$\mathit{\Omega_{inh}}$ from $i=1$ to $m$. The inherent strain component is
determined by a high-fidelity analysis such as a thermo-mechanical analysis or
by partial cutting a cantilever-shaped specimen and measuring the resulting
deformation. For details on the algorithm of the building process and
identification of inherent strain $\bm{\varepsilon}^{inh}$ using FEM, please
refer to Miki et al. [12]. The residual stress $\bm{\sigma}$ and deformation
$\bm{u}$ after completing the building process are obtained by summing the
layer-wise residual stress $\bm{\sigma}_{i}$ and deformation $\bm{u}_{i}$ as
follows:

$\displaystyle\bm{u}$ $\displaystyle=\sum_{i=1}^{m}\bm{u}_{i}\hskip
14.22636pt\text{ for }\bm{x}\in\mathit{\Omega_{p}},$ (10)
$\displaystyle\hat{\bm{\varepsilon}}^{el}$
$\displaystyle=\sum_{i=1}^{m}(\bm{\varepsilon}(\bm{u_{i}})-\bm{\varepsilon}^{inh})\hskip
14.22636pt\text{ for }\bm{x}\in\mathit{\Omega_{p}},$ (11)
$\displaystyle\bm{\sigma}$
$\displaystyle=\sum_{i=1}^{m}\bm{\sigma}_{i}=\mathbb{D}\hat{\bm{\varepsilon}}^{el}\hskip
14.22636pt\text{ for }\bm{x}\in\mathit{\Omega_{p}}.$ (12)

Based on these equations, the inherent strain $\hat{\bm{\varepsilon}}^{inh}$
on the entire part $\mathit{\Omega}$ can be derived as follows:

$\hat{\bm{\varepsilon}}^{inh}=\bm{\varepsilon}(\bm{u})-\hat{\bm{\varepsilon}}^{el}.$
(13)

Using the inherent strain $\hat{\bm{\varepsilon}^{inh}}$, the displacement
field $\hat{\bm{u}}$ of the final state is expressed by the following
governing equation:

$\displaystyle\begin{cases}\hskip 5.69054pt-\mbox{\rm
div}(\bm{\sigma})=0\hskip 14.22636pt&\text{in}\hskip
5.69054pt\mathit{\Omega}\\\ \hskip
5.69054pt\bm{\sigma}=\mathbb{D}(\bm{\varepsilon}(\hat{\bm{u}})-\hat{\bm{\varepsilon}}^{inh})\hskip
14.22636pt&\\\ \hskip 5.69054pt\hat{\bm{u}}=\bm{0}\hskip
14.22636pt&\text{on}\hskip 5.69054pt\mathit{\Gamma_{u}}\\\ \hskip
5.69054pt-\bm{\sigma}\cdot\bm{n}=\bm{0}\hskip 14.22636pt&\text{on}\hskip
5.69054pt\partial\mathit{\Omega}\setminus\mathit{\Gamma_{u}},\end{cases}$ (14)

For deriving the design sensitivity, previous studies solved the adjoint
equations corresponding to the governing equations of each intermediate state,
i.e., Eq. 8, whereas this study solves only the adjoint problem corresponding
to the governing equation of the final state, i.e., Eq. 14. In other words,
the calculation for the number of layers is replaced by a single calculation.
The deformation during the building process must be reduced to manufacture
parts with high precision in AM. In this study, we use the P-norm function to
define the deformation constraint as follows:

$G_{u}=\left(\int_{\mathit{\Omega}}\left|\hat{\bm{u}}\right|^{b}\mathrm{~{}d}\mathit{\Omega}\right)^{1/b},\\\
$ (15)

where $b\geq 2$ represents the penalization parameter set to $5$.

The design sensitivity of the AM deformation constraint $G_{u}$ is performed
based on the adjoint variable method. Denoting the adjoint variable
corresponding to the state variable $\hat{\bm{u}}$ by $\bm{\lambda}\in
H^{1}(\mathit{\Omega})^{N}$, the adjoint equation for $G_{u}$ defined as
follows:

$\displaystyle\begin{cases}\hskip 5.69054pt-\mbox{\rm
div}(\mathbb{D}\bm{\varepsilon}(\bm{\lambda}))=-\left(\int_{\mathit{\Omega}}\left|\hat{\bm{u}}\right|^{b}\mathrm{~{}d}\mathit{\Omega}\right)^{1/b-1}\left|\hat{\bm{u}}\right|^{b-2}\hat{\bm{u}}\hskip
14.22636pt&\text{in}\hskip 5.69054pt\mathit{\Omega},\\\ \hskip
5.69054pt\bm{\lambda}=\bm{0}\hskip 14.22636pt&\text{on}\hskip
5.69054pt\mathit{\Gamma_{u}},\\\ \hskip
5.69054pt-(\mathbb{D}\bm{\varepsilon}(\bm{\lambda}))\cdot\bm{n}=\bm{0}\hskip
14.22636pt&\text{on}\hskip
5.69054pt\partial\mathit{\Omega}\setminus\mathit{\Gamma_{u}},\end{cases}$ (16)

Then, the design sensitivity is derived as follows:

$G_{u}^{\prime}=-\bm{\varepsilon}(\hat{\bm{u}}):\mathbb{D}:\bm{\varepsilon}(\bm{\lambda})+\bm{\varepsilon}^{inh}:\mathbb{D}:\bm{\varepsilon}(\bm{\lambda}).$
(17)

In addition, for stabilizing the optimization calculation, the sensitivity is
averaged for each iteration process as follows:

$\tilde{G}_{u}^{\prime}=(1-\mu)\tilde{G}_{u}^{\prime old}+\mu G_{u}^{\prime},$
(18)

where $\mu$ represents a parameter that contributes to convergence and
$\tilde{G}_{u}^{\prime old}$ represents the average sensitivity of the
previous iteration. The parameter $\mu$ is set to $0.9$.

##  4 Self-support constraint

This chapter describes self-support constraint in a level-set based framework.
This self-support constraint consists of an overhang angle constraint and a
downward convex shape constraint.

###  4.1 Overhang angle constraint

In this study, a Helmholtz-type PDE filter is introduced for smoothing the
characteristic function $\chi_{\phi}$ at the structural boundary. The state
variable $\psi\in H^{1}(D)$ and its governing equation are defined as follows:

$\displaystyle\begin{cases}-aL^{2}\nabla^{2}\psi+\psi=\chi_{\phi}\hskip
14.22636pt&\text{in}\hskip 8.5359ptD\\\ \bm{n}\cdot\nabla\psi=0\hskip
109.54327pt&\text{on}\hskip 8.5359pt\partial\mathit{\Omega},\end{cases}$ (19)

where $a\in\mathbb{R}^{+}$, $L$ and $\bm{n}$ represent the diffusion
coefficient that affects the transition width of $\psi$ at the structural
boundary, representative length, and outward normal vector, respectively. The
coefficient $a$ is set to $1\times 10^{-4}$. Then, the overhang angle $\theta$
is evaluated by the inner product of the normal vector of the structural
boundary $\partial\mathit{\Omega}$ and building direction vector $\bm{d}$, as
shown in Fig. 2.

Figure 2: Conceptual diagram of the overhang angle evaluation.

The overhanging area below the minimum overhang angle $\theta_{0}$ is
represented as follows:

$\displaystyle\nabla\psi\cdot\bm{d}\geq|\nabla\psi|\cos\theta_{0}.$ (20)

From Eq. 20, the condition for constraining the overhang angle is given by

$\displaystyle{G_{o}=\int_{D}\sqrt{a}L\mathscr{R}(\nabla\psi\cdot\bm{d}-|\nabla\psi|\cos\theta_{0})\mathrm{d}\Omega},$
(21)

where $\mathscr{R}(s)$ represents the ramp function defined as follows:

$\mathscr{R}(s)=\begin{cases}0&\text{for }s\leq\epsilon_{r}\\\ s&\text{for
}s>\epsilon_{r}\end{cases}$ (22)

where $\epsilon_{r}$ represents a small parameter that limits the area below
the minimum overhang angle to only the neighborhood of the structural
boundary. The parameter $\epsilon_{r}$ is set to $1\times 10^{-2}$. Then,
denoting $p\in H^{1}(D)$ as the adjoint variable corresponding to the state
variable $\psi$, the adjoint equation for the overhang angle constraint
$G_{o}$ is defined as follows:

$\displaystyle\begin{cases}-aL^{2}\nabla^{2}p+p=\nabla\cdot\left[\sqrt{a}L\mathscr{H}(\psi)(\bm{d}-\frac{\nabla\psi}{|\nabla\psi|}\cos\theta_{0})\right]&\text{in}\hskip
8.5359ptD\\\ \bm{n}\cdot\nabla p=0\hskip 109.54327pt&\text{on}\hskip
8.5359pt\partial\mathit{\Omega},\end{cases}$ (23)

where $\mathscr{H}(s):=\mathrm{d}\mathscr{R}(s)/\mathrm{d}s$ represents the
Heviside function defined as follows:

$\mathscr{H}(s)=\begin{cases}0&\text{for }s\leq\epsilon_{r}\\\ 1&\text{for
}s>\epsilon_{r}\end{cases}$ (24)

The design sensitivity of $G_{o}$ is derived as follows:

$G_{o}^{\prime}=-p\chi_{\phi}.$ (25)

###  4.2 Downward convex shape constraint

The aforementioned overhang angle constraint creates downward convex shapes
that cannot be manufactured. This chapter describes how to suppress the
creation of such shapes by introducing a simple thermal model that can
simulate the building process [30]. Fig. 3 shows a thermal model that provides
the heat flux to the downward structural boundary layer-by-layer.

Figure 3: Domains and boundaries of the thermal model in the building process.

Similar to the inherent strain-based AM process model, part $\mathit{\Omega}$
is divided into $n$ layers to represent the building process. Domain
$\mathit{\Omega}_{T}$ is defined by each domain $\mathit{\Omega}_{j}$ for
$1\leq j\leq n$ and the domain $\mathit{\Omega}_{q}$ is determined by the
status of $\mathit{\Omega}_{j}$. Then, the temperature field $T_{j}\in
H^{1}(\mathit{\Omega}_{T})$ with heat flux $Q$ applied to the overhang
boundary $\mathit{\Gamma_{q}}$ of the added layer $\Omega_{q}$ is represented
by the following governing equation:

$\left\\{\begin{aligned} \hskip 5.69054pt&\operatorname{div}(k\nabla
T_{j})=0&\hskip 14.22636pt&\text{ in }&&\mathit{\Omega_{T}}\\\ \hskip
5.69054pt&(k\nabla T_{j})\cdot\bm{n}=Q&\hskip 14.22636pt&\text{ on
}&&\mathit{\Gamma_{q}}\\\ \hskip 5.69054pt&T_{j}=T_{amb}&\hskip
14.22636pt&\text{on }&&\mathit{\Gamma_{T}},\\\ \end{aligned}\right.$ (26)

for all indices $j=1,2,\ldots,n$, where $T_{amb}$ represents the base plate
temperature. In this thermal model, high temperatures result from shapes
obstructing heat flow. To mitigate this, considering the $\mathit{\Omega_{q}}$
heat dissipation of each layer can supress the creation of such shapes. The
downward convex shape constraint is defined as follows:

$G_{t}=\sum_{j=1}^{n}\int_{\mathit{\Omega}_{q}}\left(T_{j}-T_{amb}\right)^{2}\mathrm{~{}d}\mathit{\Omega}.$
(27)

Then, denoting the adjoint variable corresponding to the state variable
$T_{j}$ by $q_{j}\in H^{1}(\mathit{\Omega}_{T})$, the adjoint equation for the
downward convex constraint $G_{t}$ is defined as follows:

$\left\\{\begin{aligned} \hskip 5.69054pt&\operatorname{div}(k\nabla
q_{j})=2\left(T_{j}-T_{amb}\right)&\hskip 14.22636pt&\text{ in
}&&\mathit{\Omega_{T}},\\\ \hskip 5.69054pt&(k\nabla
q_{j})\cdot\bm{n}=0&\hskip 14.22636pt&\text{ on
}&&\mathit{\Omega_{T}}\setminus\mathit{\Gamma_{T}},\\\ \hskip
5.69054pt&q_{j}=T_{amb}&\hskip 14.22636pt&\text{on }&&\mathit{\Gamma_{T}},\\\
\end{aligned}\right.$ (28)

Subsequently, the design sensitivity of $G_{t}$ is derived as follows:

$G_{t}^{\prime}=\sum_{j=1}^{n}-k\nabla T_{j}\cdot\nabla q_{j}.$ (29)

##  5 Formulation of the optimization problem

In this section, we consider the minimum mean compliance problem with the
aforementioned constraints. The material domain $\mathit{\Omega}$ comprises of
isotropic linear elastic materials. The displacement is fixed at the boundary
$\mathit{\Gamma_{v}}$ and traction $\bm{t}$ is imposed on the boundary
$\mathit{\Gamma_{t}}$. In the equilibrium state, denoting the displacement
field as $\bm{v}\in H^{1}(\mathit{\Omega})^{N}$, the objective function is
given as follows:

$J_{v}=\int_{\mathit{\Gamma_{t}}}\bm{t}\cdot\bm{v}\mathrm{d}\mathit{\Gamma}.$
(30)

Under the conditions, the topology optimization problem is formulated as
follows:

$\displaystyle\inf_{\phi}\hskip 48.36967pt$ $\displaystyle J=J_{v}+\alpha
G_{u}+\beta G_{o}+\gamma G_{t}$ (31) $\displaystyle\text{subject to}:\hskip
5.69054pt$ $\displaystyle
E_{v}=\int_{\Gamma_{t}}\bm{t}\cdot\tilde{\bm{v}}\mathrm{~{}d}\Gamma-\int_{\mathit{\Omega}}\varepsilon(\bm{v}):\mathbb{D}:\varepsilon(\tilde{\bm{v}})\mathrm{~{}d}\mathit{\Omega}=0$
(32) $\displaystyle\hskip 14.22636pt\text{ for
}\forall\tilde{\bm{v}}\in\mathcal{V},\bm{v}\in\mathcal{V}$ $\displaystyle
E_{u}=\int_{\mathit{\Omega}}\hat{\bm{\varepsilon}}^{inh}:\mathbb{D}:\varepsilon(\tilde{\bm{u}})\mathrm{~{}d}\mathit{\Omega}-\int_{\mathit{\Omega}}\varepsilon(\hat{\bm{u}}):\mathbb{D}:\varepsilon(\tilde{\bm{u}})\mathrm{~{}d}\mathit{\Omega}=0$
(33) $\displaystyle\hskip 14.22636pt\text{ for
}\forall\tilde{\bm{u}}\in\mathcal{U},\hat{\bm{u}}\in\mathcal{U}$
$\displaystyle
E_{\psi}=-\int_{D}aL^{2}\nabla\psi\cdot\nabla\tilde{\psi}\mathrm{~{}d}\mathit{\Omega}-\int_{D}\psi\tilde{\psi}\mathrm{~{}d}\mathit{\Omega}-\int_{D}\chi_{\phi}\tilde{\psi}\mathrm{~{}d}\mathit{\Omega}=0$
(34) $\displaystyle\hskip 14.22636pt\text{ for
}\forall\tilde{\psi}\in\mathcal{S},\psi\in\mathcal{S}$ $\displaystyle
E_{t}=\int_{\Gamma_{q}}Q\tilde{T_{j}}\mathrm{~{}d}\Gamma-\int_{\mathit{\Omega}_{T}}\nabla
T_{j}\cdot\nabla\tilde{T_{j}}\mathrm{~{}d}\mathit{\Omega}=0$ (35)
$\displaystyle\hskip 14.22636pt\text{ for
}\forall\tilde{T_{j}}\in\mathcal{T},T_{j}\in\mathcal{T}$ $\displaystyle
G=\frac{\int_{D}\chi_{\phi}\hskip
2.84526pt\mathrm{~{}d}\mathit{\Omega}}{\int_{D}\hskip
2.84526pt\mathrm{~{}d}\mathit{\Omega}}-V_{\text{max}}\leq 0$

for all indices $j=1,2,\ldots,n$, where $\alpha$, $\beta$, and $\gamma$
represent the penalization parameters. $G$ and $V_{max}$ represent the volume
constraint and upper limit of the volume fraction, respectively. Furthermore,
functional spaces $\mathcal{V}$, $\mathcal{U}$, $\mathcal{S}$, and
$\mathcal{T}$ are defined as follows:

$\displaystyle\mathcal{V}:=\left\\{\tilde{\bm{v}}\in
H^{1}(\mathit{\Omega})^{N},\hskip 2.84526pt\tilde{\bm{v}}=\bm{0}\hskip
2.84526pt\text{ on }\hskip 2.84526pt\mathit{\Gamma}_{v}\right\\}$ (36)
$\displaystyle\mathcal{U}:=\left\\{\tilde{\bm{u}}\in
H^{1}(\mathit{\Omega})^{N},\hskip 2.84526pt\tilde{\bm{u}}=\bm{0}\hskip
2.84526pt\text{ on }\hskip 2.84526pt\mathit{\Gamma}_{u}\right\\}$ (37)
$\displaystyle\mathcal{S}:=\left\\{\tilde{\psi}\in H^{1}(\mathit{D})\right\\}$
(38) $\displaystyle\mathcal{T}:=\left\\{\tilde{T_{j}}\in
H^{1}(\mathit{\Omega_{T}}),\hskip 2.84526pt\tilde{T_{j}}=T_{amb}\hskip
2.84526pt\text{ on }\hskip 2.84526pt\mathit{\Gamma}_{T}\right\\}$ (39)

Since the minimum mean compliance problem is a self-adjoint problem, the
design sensitivity is given as follows:

$J_{v}^{\prime}=-\bm{\varepsilon}(\bm{v}):\mathbb{D}:\bm{\varepsilon}(\bm{v}).\\\
$ (40)

##  6 Numerical implementation

###  6.1 Optimization algorithm

The optimization algorithm is given below.

Step 1.

The level-set function $\phi$ is initialized.

Step 2.

The governing equations for each state variable are solved by the FEM. In the
AM process model, after solving Eq. 8, the inherent strain of the final state
is calculated and Eq. 14 is solved.

Step 3.

The objective function is computed, and if it converges, the optimization
procedure ends.

Step 4.

The adjoint equation for each adjoint variable defined is solved using the
FEM, and the design sensitivity is computed.

Step 5.

The level-set function is updated using the reaction–diffusion equation
defined in Eq. 5, and then, the procedure returns to Step 2.

The open-source PDE solver FreeFEM++[35] is used to solve the governing
equations.

###  6.2 Numerical scheme for the governing equation

We describe the ersatz material approach [36] for solving the governing
equations. In this study, a small-scaled material property is used to
represent the void domain. The boundary between the material and void domains
can be smoothly transitioned using a Heaviside function. Subsequently, the
extended elastic tensor $\tilde{\mathbb{D}}$ and thermal conductivity
$\tilde{k}$ in the fixed design domain $D$ is redefined as follows:

$\displaystyle\tilde{\mathbb{D}}(\phi;w)=\left\\{(1-c)H_{\phi}(\phi;w)+c\right\\}\mathbb{D}$
(41)
$\displaystyle\tilde{k}(\phi;w)=\left\\{(1-c)H_{\phi}(\phi;w)+c\right\\}k,$
(42)

where $H_{\phi}(\phi;w)$ is defined as

$H_{\phi}(\phi;w):=\left\\{\begin{array}[]{ll}1&\text{ for }\phi>w,\\\
\frac{1}{2}+\frac{\phi}{w}\left(\frac{15}{16}-\frac{\phi^{2}}{w^{2}}\left(\frac{5}{8}-\frac{3}{16}\frac{\phi^{2}}{w^{2}}\right)\right)&\text{
for }-w\leq\phi\leq w,\\\ 0&\text{ for }\phi<-w,\end{array}\right.$ (43)

where $w$ represents the width of the transition and $c$ represents the
scaling parameter of the material properties. Furthermore, Eq. 34 is solved by
replacing the characteristic function with this Heaviside function. In our
implementation, $c$ is set to $1.0\times 10^{-3}$, $w$ in Eqs. 32, 33, and 35
is set to $0.5$, and $w$ in Eqs. 34 is set to $0.9$.

Next, an approximate solution method for the heat conduction equation defined
in Eq. 35 is introduced. The boundary integral is replaced the domain integral
as follows:

$\int_{\mathit{\Gamma}}\xi(\bm{x})\mathrm{d}\mathit{\Gamma}\approx\int_{\mathit{\Omega}}\xi(\bm{x})\frac{\mathrm{d}H_{\psi}(\psi;w)}{\mathrm{d}\psi}|\nabla\psi|\mathrm{d}\mathit{\Omega}.$
(44)

Then, using the inner product of the normal vector
$\bm{n}_{\psi}:=\frac{\nabla\psi}{|\nabla\psi|}$ and building direction
$\bm{d}$, the boundary integral is applied to only the overhang boundary
$\mathit{\Gamma}_{q}$ as follows.

$\int_{\mathit{\Gamma}_{q}}\xi(\bm{x})\mathrm{d}\mathit{\Gamma}\approx\int_{\mathit{\Omega}_{q}}\xi(\bm{x})\frac{\mathrm{d}H_{\psi}(\psi;w)}{\mathrm{d}\psi}|\nabla\psi|\mathscr{H}(\bm{n}_{\psi}\cdot\bm{d})\mathrm{d}\mathit{\Omega}.$
(45)

##  7 Numerical examples

###  7.1 AM deformation constraint

In this section, we verify the effectiveness of the proposed method for
reducing the computational cost in the AM deformation constraint. We compare
the previous approach of solving the adjoint equations for each intermediate
state of the AM building process with the proposed approach that focuses on
only the final state. In all numerical examples, a Young’s modulus, Poisson’s
ratio, and thermal conductivity of 75 GPa, 0.34, 119 W/mK, respectively. In
Eq. 8, the in-plane component of the inherent strain $\bm{\varepsilon}^{inh}$
is set to $-0.0025$, and the penalization parameter $\alpha$ related to the AM
deformation constraint $G_{u}$ is set to $0.1$. Here, we consider a two-
dimensional cantilever beam shown in Fig. 4. The upper limit of the volume
fraction $V_{max}$ is set to $0.5$, and the representative length $L$ is set
to $25$. In addition, the fixed design domain $D$ is divided in the building
direction into $25$ layers. The model is discretized by $10,000$ elements with
an element size of $1\times 1$ mm. All numerical examples were executed on a
workstation equipped with Intel Xeon E5-2687W cores and $512$ GB of RAM.

Figure 4: Problem setting for the two-dimensional cantilever beam.

(a) Without constraint; $J_{v~{}ref}=2.52$ (b) With AM deformation constraint
(Previous method); $J_{v}/J_{v~{}ref}=103\%$

(c) With AM deformation constraint (Proposed method);
$J_{v}/J_{v~{}ref}=103\%$

Figure 5: Optimal structures with and without AM deformation constraint.
Figure 6: AM deformation for each optimal structure. Figure 7: Computational
time for each process. Figure 8: Iteration history of the objective value.

Fig. 5 illustrates optimal structures obtained through different adjoint
analysis methods. The previous method, which involves adjoint analysis for
each intermediate state of the building process, is shown in Fig. 5(b). The
proposed method, applying adjoint analysis solely to the final state, is
depicted in Fig. 5(c). Both optimal structures with the AM deformation
constraint have structures that support the overhang area and have the same
effect on the compliance. Fig. 6 shows the results of evaluating the
deformation induced by the AM building process. Figs. 6(b) and 6(c)
demonstrate that the proposed method suppresses deformation as effectively as
the previous method. Fig. 7 shows the computational time for each process in
the optimization loop. The proposed method increases the computational time
for the AM process analysis slightly; however, it significantly reduces the
computational time for the adjoint analysis. Fig. 8 shows the iterative
history of the objective value. The results indicate that the behavior of the
objective value changes slightly with a proposed method.

###  7.2 AM deformation and self-support constraints

This section presents a three-dimensional numerical example that adds the
self-support constraint to the AM deformation constraint. Fig. 9 shows the
fixed design domain and boundary conditions for a three-dimensional cantilever
beam.

Figure 9: Problem setting for the three-dimensional cantilever beam.

The upper limit of the volume fraction $V_{max}$ is set to $0.2$. The
representative length in Eq. 19 is set to $L=50$. The minimum overhang angle
is set to $\theta_{0}=45^{\circ}$, and each parameter related to Eq. 26 is set
as follows. The heat flux $Q$ is set to 10 W and the base-plate temperature
$T_{amb}$ is set to $0/,^{\circ}C$. The number of layers in the fixed design
domain is set to $m=50$ and $n=25$ for the AM deformation and downward convex
constraints, respectively. The model is discretized by $1,500,000$ elements
with an element size of $1\times 1\times 1$ mm. In addition, the penalization
parameters $\beta$ and $\gamma$ the overhand angle and downward convex
constraints are set to $0.1$ and $0.4$, respectively. The other parameters are
those set in the previous numerical example.

(a) Bird’s-eye view (b) Rear view

Figure 10: Optimal structure without constraints; $J_{v~{}ref}=1.09$.

(a) Bird’s-eye view (b) Rear view

Figure 11: Optimal structure with the self-support and AM deformation
constraints (Previous method); $J_{v}/J_{v~{}ref}=119\%$.

(a) Bird’s-eye view (b) Rear view

Figure 12: Optimal structure with the self-support and AM deformation
constraints (Proposed method); $J_{v}/J_{v~{}ref}=118\%$. Figure 13:
Comparison of the deformation induced by the building process for the
optimized 3D cantilever beam. Figure 14: Computational time for each process.
Figure 15: Iteration history of the objective value.

Figs. 10-12 illustrate the optimal structures under different conditions:
without constraints (Fig. 10), with constraints imposed by the previous method
(Fig. 11), and with constraints imposed by the proposed method (Fig. 12). The
results show that, compared to the unconstrained optimal structure, both
constrained optimal structures are self-supporting shapes with the same
geometric features. Fig. 13 presents the results of the deformations analyzed
for optimal structures using the AM process model. Figs. 13(b) and 13(c)
illustrate that the proposed method suppresses deformation effectively,
comparable to the performance of the previous method. Fig. 14 shows the
computational time for each process in the optimization loop. Compared to the
two-dimensional numerical result, the computational time for the adjoint
analysis is reduced significantly. Because the number of elements is relevant.
Therefore, the larger the model with a large number of elements, the more
effectively the proposed method can reduce the computational time. Fig. 15
shows the iterative history of the objective function, which indicates that
the convergence behavior remains unchanged even with the inclusion of the
self-support constraint. Thus, given these results, the proposed method can
effectively reduce the computational time of the AM deformation constraints.

##  8 Conclusion

This paper proposed a computational cost reduction method for the AM
deformation constraint in topology optimization. The main contributions of
this study are summarized as follows:

  1. 1.

An AM process model based on the inherent strain method was presented, and the
governing equation for the final state of the building process were
formulated.

  2. 2.

A constraint function that suppressed the AM deformation was formulated. Then,
we proposed a method to reduce the computational cost of sensitivity analysis
using the governing equations that represent the final state of the building
process.

  3. 3.

We formulated an optimization problem with the AM deformation constraint
combined with the self-support constraint and constructed an algorithm using
FEM.

  4. 4.

Through numerical examples, we confirmed that the proposed method reduced the
computation time for optimization and obtained an optimal structure with the
same performance as the previous method.

## Acknowledgments

This study was supported by a JSPS grant for Scientific Research (C)
JP21K03826.

## Statement

During the preparation of this work the author used DeepL and ChatGPT in order
to English proofreading. After using this tool, the authors reviewed and
edited the content as needed and takes full responsibility for the content of
the publication.

## References

## References

  * [1] A. T. Gaynor, J. K. Guest, Topology optimization considering overhang constraints: Eliminating sacrificial support material in additive manufacturing through design, Structural and Multidisciplinary Optimization 54 (5) (2016) 1157–1172.
  * [2] M. Langelaar, Topology optimization of 3d self-supporting structures for additive manufacturing, Additive manufacturing 12 (2016) 60–70.
  * [3] M. Langelaar, An additive manufacturing filter for topology optimization of print-ready designs, Structural and multidisciplinary optimization 55 (3) (2017) 871–883.
  * [4] G. Allaire, C. Dapogny, R. Estevez, A. Faure, G. Michailidis, Structural optimization under overhang constraints imposed by additive manufacturing technologies, Journal of Computational Physics 351 (2017) 295–328.
  * [5] X. Qian, Undercut and overhang angle control in topology optimization: a density gradient based integral approach, International Journal for Numerical Methods in Engineering 111 (3) (2017) 247–272.
  * [6] Y. Wang, J. Gao, Z. Kang, Level set-based topology optimization with overhang constraint: towards support-free additive manufacturing, Computer Methods in Applied Mechanics and Engineering 339 (2018) 591–614.
  * [7] E. van de Ven, R. Maas, C. Ayas, M. Langelaar, F. van Keulen, Continuous front propagation-based overhang control for topology optimization with additive manufacturing, Structural and Multidisciplinary Optimization 57 (5) (2018) 2075–2091.
  * [8] A. Garaigordobil, R. Ansola, J. Santamaría, I. Fernández de Bustos, A new overhang constraint for topology optimization of self-supporting structures in additive manufacturing, Structural and Multidisciplinary Optimization 58 (5) (2018) 2003–2017.
  * [9] C.-J. Thore, H. A. Grundström, B. Torstenfelt, A. Klarbring, Penalty regulation of overhang in topology optimization for additive manufacturing, Structural and Multidisciplinary Optimization 60 (1) (2019) 59–67.
  * [10] R. A. Wildman, A. T. Gaynor, Topology optimization for reducing additive manufacturing processing distortions, Tech. rep., Weapons and Materials Research Directorate, US Army Research Laboratory … (2017).
  * [11] G. Allaire, L. Jakabčin, Taking into account thermal residual stresses in topology optimization of structures built by additive manufacturing, Mathematical Models and Methods in Applied Sciences 28 (12) (2018) 2313–2366.
  * [12] T. Miki, T. Yamada, Topology optimization considering the distortion in additive manufacturing, Finite Elements in Analysis and Design 193 (2021) 103558\.
  * [13] G. Allaire, B. Bogosel, Optimizing supports for additive manufacturing, Structural and Multidisciplinary Optimization 58 (6) (2018) 2493–2515.
  * [14] M. Zhou, Y. Liu, Z. Lin, Topology optimization of thermal conductive support structures for laser additive manufacturing, Computer Methods in Applied Mechanics and Engineering 353 (2019) 24–43.
  * [15] C. Wang, X. Qian, Optimizing support for heat dissipation in additive manufacturing, in: International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, Vol. 83983, American Society of Mechanical Engineers, 2020, p. V009T09A018.
  * [16] T. Miki, S. Nishiwaki, Topology optimization of the support structure for heat dissipation in additive manufacturing, Finite Elements in Analysis and Design 203 (2022) 103708.
  * [17] R. Ranjan, C. Ayas, M. Langelaar, F. v. Keulen, Controlling local overheating in topology optimization for additive manufacturing, Structural and Multidisciplinary Optimization 65 (6) (2022) 1–21.
  * [18] L. Papadakis, A. Loizou, J. Risse, J. Schrage, Numerical computation of component shape distortion manufactured by selective laser melting, Procedia Cirp 18 (2014) 90–95.
  * [19] N. Hodge, R. Ferencz, J. Solberg, Implementation of a thermomechanical model for the simulation of selective laser melting, Computational Mechanics 54 (1) (2014) 33–51.
  * [20] C. Li, C. Fu, Y. Guo, F. Fang, A multiscale modeling approach for fast prediction of part distortion in selective laser melting, Journal of Materials Processing Technology 229 (2016) 703–712.
  * [21] T. Mukherjee, W. Zhang, T. DebRoy, An improved prediction of residual stresses and distortion in additive manufacturing, Computational Materials Science 126 (2017) 360–372.
  * [22] E. R. Denlinger, M. Gouge, J. Irwin, P. Michaleris, Thermomechanical model development and in situ experimental validation of the laser powder-bed fusion process, Additive Manufacturing 16 (2017) 73–80.
  * [23] M. Chiumenti, E. Neiva, E. Salsi, M. Cervera, S. Badia, J. Moya, Z. Chen, C. Lee, C. Davies, Numerical modelling and experimental validation in selective laser melting, Additive Manufacturing 18 (2017) 171–185.
  * [24] C. Li, J. Liu, X. Fang, Y. Guo, Efficient predictive model of part distortion and residual stress in selective laser melting, Additive Manufacturing 17 (2017) 157–168.
  * [25] N. Keller, V. Ploshikhin, New method for fast predictions of residual stress and distortion of am parts, in: Solid Freeform Fabrication Symposium, Vol. 25, 2014.
  * [26] I. Setien, M. Chiumenti, S. van der Veen, M. San Sebastian, F. Garciandía, A. Echeverría, Empirical methodology to determine inherent strains in additive manufacturing, Computers & Mathematics with Applications 78 (7) (2019) 2282–2295.
  * [27] Q. Chen, X. Liang, D. Hayduke, J. Liu, L. Cheng, J. Oskin, R. Whitmore, A. C. To, An inherent strain based multiscale modeling framework for simulating part-scale residual deformation for direct metal laser sintering, Additive Manufacturing 28 (2019) 406–418.
  * [28] X. Liang, Q. Chen, L. Cheng, D. Hayduke, A. C. To, Modified inherent strain method for efficient prediction of residual deformation in direct metal laser sintered components, Computational Mechanics 64 (6) (2019) 1719–1733.
  * [29] B. C. Prabhune, K. Suresh, A fast matrix-free elasto-plastic solver for predicting residual stresses in additive manufacturing, Computer-Aided Design (2020) 102829.
  * [30] T. Miki, Self-support topology optimization considering distortion for metal additive manufacturing, Computer Methods in Applied Mechanics and Engineering 404 (2023) 115821.
  * [31] S. Xu, J. Liu, Y. Ma, Residual stress constrained self-support topology optimization for metal additive manufacturing, Computer Methods in Applied Mechanics and Engineering 389 (2022) 114380.
  * [32] T. Yamada, K. Izui, S. Nishiwaki, A. Takezawa, A topology optimization method based on the level set method incorporating a fictitious interface energy, Computer Methods in Applied Mechanics and Engineering 199 (45-48) (2010) 2876–2891.
  * [33] M. Bihr, G. Allaire, X. Betbeder-Lauque, B. Bogosel, F. Bordeu, J. Querois, Part and supports optimization in metal powder bed additive manufacturing using simplified process simulation, Computer Methods in Applied Mechanics and Engineering 395 (2022) 114975.
  * [34] M. Bugatti, Q. Semeraro, Limitations of the inherent strain method in simulating powder bed fusion processes, Additive Manufacturing 23 (2018) 329–346.
  * [35] F. Hecht, New development in freefem++, J. Numer. Math. 20 (3-4) (2012) 251–265.
URL https://freefem.org/

  * [36] G. Allaire, F. Jouve, A.-M. Toader, Structural optimization using sensitivity analysis and a level-set method, Journal of Computational Physics 194 (1) (2004) 363–393.


# Non-central panorama indoor dataset

Bruno Berenguel-Baeta
Instituto de Investigacion en Ingenieria de Aragon
Department of Computer Science and Systems Engineering
University of Zaragoza, Zaragoza, Spain
<EMAIL_ADDRESS>
& Jesus Bermudez-Cameo
Instituto de Investigacion en Ingenieria de Aragon
Department of Computer Science and Systems Engineering
University of Zaragoza, Zaragoza, Spain
<EMAIL_ADDRESS>
& Jose J. Guerrero
Instituto de Investigacion en Ingenieria de Aragon
Department of Computer Science and Systems Engineering
University of Zaragoza, Zaragoza, Spain
<EMAIL_ADDRESS>
Corresponding author.

###### Abstract

Omnidirectional images are one of the main sources of information for learning
based scene understanding algorithms. However, annotated datasets of
omnidirectional images cannot keep the pace of these learning based algorithms
development. Among the different panoramas and in contrast to standard central
ones, non-central panoramas provide geometrical information in the distortion
of the image from which we can retrieve 3D information of the environment [2].
However, due to the lack of commercial non-central devices, up until now there
was no dataset of these kinds of panoramas. In this data paper, we present the
first dataset of non-central panoramas for indoor scene understanding. The
dataset is composed by 2574 RGB non-central panoramas taken in around 650
different rooms. Each panorama has associated a depth map and annotations to
obtain the layout of the room from the image as a structural edge map, list of
corners in the image, the 3D corners of the room and the camera pose. The
images are taken from photorealistic virtual environments and pixel-wise
automatically annotated.

_Keywords_ Computer Vision, Indoor Scene Understanding, Non-central Panoramas,
Omnidirectional Vision, Monocular Depth Estimation, Layout Estimation

††A final version of this article can be found at
https://doi.org/10.1016/j.dib.2022.108375

Specifications Table

Subject |  Computer Science: Computer Vision and Pattern Recognition
---|---
Specific subject area |  Non-central circular panoramas for indoor scene understanding.
Type of data |  RGB Image (.png)
Color code depth maps (.png)
Layout annotations (.png, .npy, .txt, .mat)
How the data were acquired |  Random generation of virtual environments.
Ad-hoc programmable camera projection model for image rendering via ray
tracing. The RGB images are rendered with POV-Ray111The Persistence of Vision
Raytracer.http://www.povray.org (accessed May 2022) and the depth maps with
Mega-POV222MegaPOV. http://megapov.inertart.net (accessed May 2022).
Layout annotations are obtained from the 3D model of the virtual environment.
Data format |  Raw
Filtered
Description of data collection |  From the generated virtual environments, we randomly place the non-central acquisition system in different locations inside the environment. The radius of acquisition is 1 meter. The size of the panoramas is 1024x512 pixels. Once acquired the non-central panoramas, we exclude those that will be physically imposible to acquire in a real situation (i.e. the non-central acquisition system goes through an object or a wall, creating a black hole in the image).
Data source location |  Institution: University of Zaragoza, Department of Computer Science and Systems Engineering
City/Town/Region: Zaragoza, Aragon
Country: Spain
Latitude and longitude for collected samples/data: 41.68390649378424 N
-0.8887470805962938 E
Data accessibility |  Repository name: Google Drive (via Github repository to make it easier for the user)
URL to Github repository:
https://github.com/jesusbermudezcameo/NonCentralIndoorDataset
URL to data (Google Drive shared folder):
https://drive.google.com/drive/folders/18OQXpbZsr3RBphU0kJC0OS2OXr-3BrkV?usp=sharing
Related research
article |  B. Berenguel-Baeta, J. Bermudez-Cameo and J.J. Guerrero, Atlanta Scaled Layouts from Non-central Panoramas. Pattern Recognition (2022). DOI:https://doi.org/10.1016/j.patcog.2022.108740

Value of the Data

  * $\bullet$

The presented dataset is the first existing dataset with non-central
panoramas. Besides it includes annotations for different purposes as layout
recovery, line extraction and depth estimation.

  * $\bullet$

Researchers who want to take advantage of the geometrical properties of non-
central systems can find in this dataset a perfect source of information for
evaluation and development of new algorithms.

  * $\bullet$

Since it is the first non-central dataset, it can be used to adapt existing
algorithms for omnidirectional central images to the no-central case. Besides,
in the related research [1], only the RGB images and layout annotations have
been used, leaving the depth maps for future research topics.

Data Description
The dataset contains a set of gravity oriented panoramas. We make this
specification since non-central panoramas cannot be rotated as a data
augmentation in a different axis that the revolution axis of the non-central
system. The dataset includes the folders: img contains the RGB non-central
panoramas are located; depth_coded contains depth maps coded in 3 channels
(RGB channels) associated with the RGB panoramas; EM_gt contains one channel
images where the structural lines of the environments are defined; DataPython
and DataMatLab contains ground truth information used to evaluate the work
presented in [1], including the 3D position of the corners of the room
(3D_gt), the camera location (cam_pose), the labelling for the floor-wall and
ceiling-wall intersections as spherical coordinates of the projecting rays
(label_ang) and the pixel coordinates of the 3D corners in the panoramas
(label_cor).

The main folder follows the following distribution.

  * •

NonCentralIndoorDataset

    * –

img

    * –

depth_coded

    * –

EM_gt

    * –

DataPython

      * *

3D_gt

      * *

cam_pose

      * *

label_ang

      * *

label_cor

    * –

DataMatLab

      * *

mat_gt

Experimental design, materials and methods

The data of this dataset is obtained from environments randomly and
synthetically generated. We first generate a random layout constrained by
different structural limits (see Table 2). These limits include minimum and
maximum: area of the room, number of walls, wall length, angle between walls
(if non-Manhattan), room height and Manhattan ratio.

Table 2: Layout limits for randomized generation. Parameter | Min. | Max.
---|---|---
Area $[m^{2}]$ | 25 | 110
n. walls | 4 | 14
wall length $[m]$ | 0.5 | 8
walls angle $[degs]$ | 25 | 100
Room height $[m]$ | 2.5 | 4.25
Manhattan ratio | 0.7

The first version of each layout is a Manhattan layout constrained by the
previous limits. Then, depending on the Manhattan ratio, a set of random
vertical planes clip the layout introducing oblique walls in angles with
respect to the previous walls in the range of the walls angle limit. After the
layout clip, we evaluate the final layout, checking if it meets the
constraints. If any of the constraints is not satisfied, the layout is deleted
and a new one is generated.

Once the structure of the room is obtained, we randomly set different kinds of
walls containing doors, windows, colors and textures. These characteristics of
the walls are randomly selected from different pools (e.g. we have different
models of doors and windows to select). At this point we select the color and
textures of the ceiling and floor of the room. Once defined the 3D structure
of the room, we include objects in it. For that purpose, we build a free space
map where we can place different kinds of objects. First we consider the kind
of objects that are placed next to a wall in a fixed orientation (beds,
wardrobes, desks). Then we place objects in the room in a random position and
orientation (chairs, sofas, carpets). Finally, we place objects that are
placed on top of other objects (cups, clocks, clothes) and we also place
lighting sources for a more realistic rendering. All of these objects are
randomly picked from different pools depending on the object class. Besides,
the ambient illumination conditions are also randomly picked from 3 different
configurations.

(a)  (b) (a)

(c)

(d)  (e) (b)

(f)

Figure 1: Layout generation and non-central panorama rendering examples. (a)
Manhattan random room. (b) Atlanta random room.

Once defined our room, we have all the necesary information to generate the
ground truth and the labelling of the dataset (see Fig. 1). The next step is
to render the RGB non-central panoramas and generate the depth maps. The color
images are rendered with the ray tracing software POV-Ray while the depth maps
are obtained with the use of MegaPOV. The non-central panoramic camera is
modeled by using an ad-hoc programmable camera projection model included in
last versions of POV-Ray. For each scene we can render different acquisitions
modifying the position and the orientation of the camera. In the case of the
proposed dataset, the camera is always oriented with the gravity direction.
Notice that, by contrast with central panoramas (e.g. equirectangular images),
we can not post-process a single render for obtaining different panoramas with
different orientations.

Acknowledgments
Funding: This work was supported by RTI2018-096903-B-100 (AEI/ FEDER, UE).

References

  * ${[1]}$

B.Berenguel-Baeta, J.Bermudez-Cameo and J.J.Guerrero. Atlanta scaled layouts
from non-central panoramas. Pattern Recognition. (2022).
https://doi.org/10.1016/j.patcog.2022.108740

  * ${[2]}$

J. Bermudez-Cameo, O. Saurer, G. Lopez-Nicolas, J.J. Guerrero and M.
Pollefeys. Exploiting line metric reconstruction from non-central circular
panoramas. Pattern Recognition Letters, 94, 30-37. (2017).
https://doi.org/10.1016/j.patrec.2017.05.006