TITLE: Find total number of balls with conditional probability
QUESTION [0 upvotes]: So, this is the problem
You have an urn with balls of two different colors and you take two without replacement. The probability of getting two of the same color is the same as getting one from each color. How many balls from each color if we know that in total we have 31?
What I do is (X,Y denote the different colors, index denotes different drawings)
$P(X_2|X_1) = P(Y_2|X_1) =
P(X_2\cap X_1)/P(X_1) = P(Y_2\cap X_1)/P(X_1)  \\ P(X_2\cap X_1) = P(Y_2\cap X_1) \rightarrow\\
\frac{X-1}{29}\frac{X}{30} = \frac{Y}{29}\frac{X}{30} \rightarrow\\
X-1=Y \rightarrow X= 16, Y=15 $
Is this correct ? I am worried, that one could choose instead of $P(Y_2\cap X_1)/P(X_1)$ to have $P(X_2\cap Y_1)/P(Y_1)$ and then maybe my solution would not hold...

REPLY [1 votes]: Two of the same colour is the event: $(X_1\cap X_2)\cup(Y_1\cap Y_2)$.
One from each colour is the event: $(X_1\cap Y_2)\cup(Y_1\cap X_2)$
Remember that there are $x$ balls of one colour and $31-x$ of the other.

REPLY [0 votes]: If there is $x$ of one color, there is $(31-x)$ of other.
Number of ways of picking balls of the same color,
$ \displaystyle x(x-1)  + (31-x) (30-x)$
Number of ways of picking one ball of each color,
$= 2 x (31-x)$
If you equate the two, you get $ \displaystyle  x = \frac{31 \pm \sqrt {31}}{2}$.
So there is in fact no way to have $31$ balls in total of two different colors such that the probability of picking two balls of same color is same as probability of picking one ball of each color (without replacement).


text:
\section{Numerical examples}\label{sec:numerical_examples}

\subsection{Implementation}
The MultiMesh FEM is implemented in the finite element framework FEniCS~\cite{alnaes2015fenics,logg2010dolfin}, and will be released in version 2018.1.0. All meshes used in this paper have been generated with GMSH, version 3.0.6~\cite{gmsh}. The complete code of the examples are published at \href{https://bitbucket.org/dokken/multimeshshapeopt}{Bitbucket} and access can be granted by emailing the corresponding author.

\subsection{Optimization of Current Carrying Multi-cables}\label{sec:multicable}

In modern cars, the number of electronic devices have been increasing,
especially in electric and hybrid cars. This means that car manufacturers must design
wires carrying current to the different devices as compactly as possible. An example of such a multi-cable is shown in \cref{fig:multicable_example}.
 This motivates optimizing the design of such multi-cable to minimize the heat inside, see~\cite{harbrecht2016} and the references therein.
As design variables, the authors chose the position of each internal cable of the multi-cable.
Each optimization iteration results in new cable positions and a re-meshing strategy was used to
update the mesh. Since the internal cables are translated within the multi-cable, the re-meshing step can be avoided if we apply \cref{alg:shape}.

\begin{figure}[!ht]
\centering
\includegraphics[width=0.2\linewidth]{multicable_tidal_turbine_nova_scotia_seen_by_funke.jpg}
\caption{A current carrying multi-cable as studied in \cref{sec:multicable}.}
\label{fig:multicable_example}
\end{figure}

To demonstrate this, we consider a simplified multi-cable optimization problem:
\begin{align}  \label{eq:multicable:functional}
  \min_{\Omega, T}J(\Omega,T)=\int_{\Omega}\frac{1}{q}\vert T\vert^q \md x,
  \quad q>1,
\end{align}
subject to
\begin{align}
  \begin{alignedat}{2}
    - \nabla \cdot (\lambda \nabla T) - c T &= f  && \text{ in } \Omega, \\
    \lambda \frac{\partial T}{\partial n} + (T-T^{\text{ex}}) &= 0 &&
    \text{ on } \partial\Omega,\\
         \llbracket T \rrbracket_{\Gamma_{i/e}^j} &= 0 && \text{ on } \Gamma_{i}^{j}\cup\Gamma_{e}^{j}, j=1,\dots, N,\\
           \left\llbracket\lambda\der{T}{n}\right\rrbracket_{\Gamma_{i/e}^j} &= 0 && \text{ on } \Gamma_{i}^{j}\cup\Gamma_{e}^{j}, j=1,\dots, N,
    \label{eq:multicable:state_system}
  \end{alignedat}
\end{align}
where $\Omega = \Omega_{fill} \cup \Omega_{insulation} \cup \Omega_{metal}$ describes a 2D slice through the multi-cable with $N$ internal cables, as specified in \cref{fig:MultiCableDomain} a). The operation $\llbracket \cdot \rrbracket_{\Gamma_{i/e}^j}$denotes the jump over the interface $\Gamma^j_i$ and $\Gamma^j_e$.
The state equation is a Poisson equation with an additional linear source term with temperature coefficient $c=0.04$. This term describes the rise of electrical resistivity for increasing temperatures in conductive material. The external boundary condition is a Robin-condition, related to the air surrounding the cable, with temperature $T^{\text{ex}}=3.2$.
Furthermore, we set $q=3$ to approximate the $L^\infty$ norm, as done by ~\cite{harbrecht2016}.
A detailed derivation of these equations can be found in~\cite{loos2014joule}.
The source-term $f$ and heat-conductivity $\lambda$ are discontinuous, piece-wise constant functions with the following values:

\begin{center}
  \footnotesize
  \begin{tabular}{|c|c|c|c|}
    \hline
      & $\Omega_{fill}$ &  $\Omega_{insulation}$ & $\Omega_{metal}$ \\
      \hline
      $f$ & 0.0 & 0.0 & 50.0 \\ \hline
      $\lambda$ & 0.08 &  0.19 & 40.0\\ \hline
\end{tabular}
\end{center}
\noindent

The state equation is solved with the MultiMesh FEM for arbitrarily
many intersecting meshes, see \cref{sec:MMFEM} and~\cite{MMFEMarbit}.
There are different ways of creating the overlapping meshes for this
problem. We chose to represent the domain by one mesh for the filling material, and $N$ meshes for the inner cables,
see \cref{fig:MultiCableDomain}b).
The meshes for the internal cables include a halo of filling material, which was chosen sufficiently large so that the heat conductivity $\lambda$ is constant over the cells categorized as overlapped. This guarantees that the solution on the overlap area does not add non-physical contributions to the weak formulation.

We choose the position of each inner cable as a design variable, as described in \cref{sec:updatemesh}.
Since the background cable is fixed, there are additional constraints on the centroids.
The distance of each cable to origin has to be bounded, such that the cables do not move
outside the cable defined by $\Gamma^{ex}$. To enforce these
constraints, we used a projected Armijo Rule~\cite{hinze2008optimization}.
\begin{figure}[!ht]
  \centering
  \begin{tikzpicture}
    \node [above right, inner sep=0pt]{
      \includegraphics[width=\linewidth]{Cable_Illustration.pdf}};
    \node (B) at (0.15\linewidth,-0.25) {(a)};
    \node (B) at (0.6\linewidth,-0.25) {(b)};
  \end{tikzpicture}
  \caption{(a) Illustration of the material composition of a multi-cable with annotated boundaries. (b) Illustration of how an internal cable is represented by a separate domain. Every domain includes an extra halo surrounding the cable. A continuity condition is therefore imposed for the state-equations at $\Lambda_1$.}
  \label{fig:MultiCableDomain}
\end{figure}

The adjoint equations are derived following \cref{sec:adjoint} and are:
\begin{align}
  \begin{alignedat}{2}
    -\nabla\cdot(\lambda\nabla p)-cp &= -T\vert T\vert^{q-2}&&
    \text{ in }\Omega,\\
    p&=p_{ex}&& \text{ on } \partial\Omega,\\
    \lambda\der{p}{n}+p_{ex}&=0 &&\text{ on }
    \partial\Omega,\\
    \left\llbracket \lambda\der{p}{n} \right\rrbracket_{\Gamma_{i/e}^j} &=0 && \text{ on } \Gamma_{i}^{j}\cup\Gamma_{e}^{j}, j=1,\dots, N,\\
    \llbracket p  \rrbracket_{\Gamma_{i/e}^j} &=0 && \text{ on } \Gamma_{i}^{j}\cup\Gamma_{e}^{j}, j=1,\dots, N.

  \end{alignedat}
\end{align}
and the shape sensitivity is
\begin{align}
  \begin{alignedat}{2}
    \md J(\Omega)[s] &= \sum_{j=1}^N \Int{\Gamma_{i}^{j}\cup\Gamma_{e}^{j}}{}s\cdot n\Big(&&
     \llbracket -cTp-fp \rrbracket_{\Gamma_{i/e}^j}
  -\lambda^+\der{p^+}{n}\left\llbracket \der{T}{n} \right\rrbracket_{\Gamma_{i/e}^j}\\
&  && +\llbracket \lambda \rrbracket_{\Gamma_{i/e}^j}\nabla_\Gamma p^+\cdot \nabla_\Gamma T^+
  \Big) \md S,
  \end{alignedat}\label{eq:MultiCableGradient}
\end{align}
where the super-script $+$ denotes the evaluation of a function from
the fill side at $\Gamma_{e}^j$, and evaluation at the insulation side of $\Gamma_{i}^j$.

\subsubsection{Results}\label{Multicable:results}
Numerical results were performed on a multi-cable with radius $1.2$, containing $N$ internal cables with $0.2$ radius not counting insulation and a $0.055$ thick insulation.

The adjoint equation and shape sensitivity were verified by a Taylor test, where one
inner cable without insulation was placed at $(0.03,0.2)$.
A MultiMesh consisting of a background mesh with $33802$ elements
and $15842$ elements for the mesh of each of the internal cable was used.
The convergence rates in the steepest descent direction are shown in \cref{tab:MultiCable1Taylorsteep}.
Similar convergence rates were obtained for other perturbations,
indicating that the adjoint equation and the shape derivatives are correct.
Further tests with different mesh refinements showed that the Taylor convergence rates
decreases if the mesh is too coarse, due to the discrete inconsistency in the functional sensitivities, see \cref{sec:shape}.

\input{taylor_MC1steep.tex}

Next, the optimization algorithm was verified by optimizing the position of three
identical cables. For this setup, it is known that the optimal positioning of the cables form an equilateral triangle~\cite{harbrecht2016}.
The optimization loop was terminated when the relative functional reduction dropped below a set tolerance, i.e.
$\left\vert\frac{J(\Omega^{k+1})-J(\Omega^k))}{J(\Omega^{k+1})}\right\vert<10^{-6}$. The optimal configuration was found after $62$ iterations, when the functional decreased from $3.2\cdot 10^4$ to $3.9\cdot 10^3$. The optimized angles between the cables were $58.7, 62.1, 59.3$ degrees. The initial and optimized configurations are shown in \cref{fig:MultiCable3}.

Furthermore, we compared the computational expense of the MultiMesh FEM with a traditional FEM approach.
For this comparison, we considered a problem with one internal cable. We measured the run-time of one approximate iteration, consisting of assembling and solving the state and adjoint systems, and updating the mesh.
For the MultiMesh-approach, the mesh updated consists of translating the mesh coordinates of the inner cable and to recompute the new mesh intersections.

For the traditional one mesh approach, the mesh update was performed through remeshing.
The linear systems arising in both approaches were solved using the an LU solver.
The timing results are shown in \cref{tab:Timing_Poisson}.
The results show that the assembly of the MultiMesh-system is more time consuming that the traditional one mesh approach, primarily caused by the the additional stabilization terms. However, this additional expense is outweighed by a significant lower cost for  the mesh update compared to re-meshing. Therefore, the estimated iteration cost for the MultiMesh-approach is a third of the traditional approach.

Finally, the optimization was performed on a problem with five internal cables of different sizes and with different insulation parameters. The parameters are listed in \cref{tab:MC5}. The initial and optimized cable configurations are shown in \cref{fig:MC5}. The results show that the smallest cable, Cable 5, is placed far away from the other cables. This can be explained through the fact that this cable has the lowest insulation parameter $\lambda_{iso}$ and the largest heat source $f$. The optimal solution was found after $103$ iterations, and the functional decreased from $1.1\cdot 10^5$ to $9.8\cdot 10^3$.

\begin{figure}[!ht]
  \centering
  \begin{tikzpicture}
    \node [above right, inner sep=0pt]{
      \includegraphics[width=\linewidth]{Cables_3.png}
    };
    \node (B) at (0.2\linewidth,-0.2) {(a)};
    \node (B) at (0.7\linewidth,-0.2) {(b)};
  \end{tikzpicture}
   \caption{
     (a) The cable configuration and temperature distribution the three cables before the optimization.
     (b) The cable configuration and temperature distribution after the optimization. The inner cables from an equilateral triangle.}\label{fig:MultiCable3}
\end{figure}

\input{Timing_Poisson.tex}

\begin{table}[!ht]
  \centering
  \caption{The setup for the 5 multi-cable optimization shown in \cref{fig:MC5}. The parameters $\lambda_{fill},\lambda_{metal}$, are the same as in the other simulations.}\label{tab:MC5}
  \footnotesize
  \begin{tabular}{|c|c|c|c|c|c|}
    \hline   & Cable 1 & Cable 2 & Cable 3 & Cable 4 & Cable 5\\
    \hline Init. Positions & $0, 0.45$ & $-0.4,-0.15$ & $0.2,-0.4$ & $0.5,0.15$ & $0.7,-0.3$  \\
    \hline Opt. Positions & $-0.434, 0.777$ & $-0.903,  0.060$  &
    $-0.191, -0.9$ &  $0.536,  0.784$  & $0.902, -0.384$  \\
    \hline $r_{iso}$ & $0.255$ & $0.242$ & $0.218$ & $0.174$ & $0.122$ \\
    \hline $r_{met}$ & $0.2$ & $0.19$ & $0.171$ & $0.137$ & $0.096$ \\
    \hline $\lambda_{iso}$  & $0.19$ & $0.162$ & $0.133$ & $0.105$ & $0.076$ \\
    \hline $f$ & $50$ & $60$ & $70$ & $80$ & $90$ \\ \hline
  \end{tabular}
\end{table}

\begin{figure}[!ht]
  \centering
  \begin{tikzpicture}
    \node [above right, inner sep=0pt]{
      \includegraphics[width=\linewidth]{Cables_5.png}
    };
    \node (B) at (0.2\linewidth,-0.2) {(a)};
    \node (B) at (0.7\linewidth,-0.2) {(b)};
  \end{tikzpicture}
  \caption{(a) Initial configuration of five cables of different sizes and with different material parameters. (b) The optimized positioning of the cables. The smallest cable is placed as far away from the other cables because it has the lowest insulation and highest heat source.} \label{fig:MC5}
\end{figure}
\FloatBarrier

\subsection{Shape Optimization of an Obstacle in Stokes Flow}\label{sec:Stokes}
This example considers the drag minimization of an object
subject to a Stokes flow in two dimensions. This problem has a known analytical solution
first presented in~\cite{pironneau1974optimum}.
The drag is measured by the dissipation of kinetic energy into heat, that is
\begin{align}
  J_S &=\Int{\Omega}{} \sum_{i,j=1}^2\left(\der{u_i}{x_j}\right)^2\md x,
\end{align}
where $\der{u_i}{x_j}$ denotes the derivative of the i-th velocity component in
the j-th direction.
The domain consists of a unit square excluding an obstacle as shown in \cref{fig:channel} (a).
\begin{figure}[ht!]
  \centering
  \begin{tikzpicture}
    \node [above right, inner sep=0pt]{
      \includegraphics[width=\linewidth]{Stokes_Illustration.pdf}
    };
    \node (B) at (0.25\linewidth,-0.3) {(a)};
    \node (B) at (0.65\linewidth,-0.3) {(b)};
  \end{tikzpicture}
  \caption{(a) The computational domain of shape optimization problem in Stokes flow with
  labeled boundaries.
  The fluid part of the domain is shown in blue, while the obstacle is shown in gray.
  (b) A coarse variant of the premeshes visualized by black and red edges. The different cell-types of the background mesh are marked in colors.
  }\label{fig:channel}
\end{figure}
The trivial solution to this problem would be to remove the object from the Stokes-flow completely.
This is avoided by introducing additional constraints on the area and
centroid of the obstacle.
Denoting the target centroid and area of the obstacle as $(C_{x0},C_{y0})=(0.5,0.5)$ and the $V_O =0.05 $ respectively,
we enforce constraints with quadratic penalty terms, yielding the cost functional
\begin{equation}
  \begin{alignedat}{1}
    J_V &= \gamma_1 \left(\vert \Omega \vert - \vert \Omega_0 \vert \right)^2,\\
    J_{Cx} &= \gamma_2 (\text{C}_x-\text{C}_{x0})^2, \\
    J_{Cy} &= \gamma_2 (\text{C}_y-\text{C}_{y0})^2,\\
    J &= J_S + J_V + J_{Cx} + J_{Cy},
    \end{alignedat}
\end{equation}
with penalty parameters $\gamma_1 > 0$ and $\gamma_2 > 0$. We denote the target fluid area as $\vert\Omega_0\vert=1-V_O=0.95$,
the actual fluid area as $\vert \Omega \vert = \int_\Omega 1 \md x$,  and
the coordinate of the obstacle's centroid as $C_x$ and $C_y$, for instance $C_x =\Big(\half-\Int{\Omega}{}x\md x \Big)
 / {\Big(1- \vert \Omega \vert \Big)}$.

To summarize, we can write the optimization problem as
\begin{align}
  \begin{split}\label{eq:Stokes:functional}
    \min_{\Omega, u} J(\Omega, u)

  \end{split}
\end{align}
subject to
\begin{equation}
  \label{eq:Stokes:state}
  \begin{alignedat}{2}
    -\Delta  u + \nabla p &= 0 &&\quad\text{ in } \Omega,\\
    \nabla\cdot  u &= 0 &&\quad\text{ in } \Omega,\\
     u &= 0 &&\quad\text{ on }\Gamma_2,\\
     u &= u_0 &&\quad\text{ on } \Gamma_1\cup\Gamma_3,\\
    \der{u}{n} + pn &= 0 &&\quad\text{ on } \Gamma_4,
  \end{alignedat}
\end{equation}
where $p$ is the fluid pressure, $u_0$ a prescribed boundary velocity, and the domain $\Omega$ is a function of $\Gamma_2$. The boundaries $\Gamma_i, i=1,\dots,4$ are visualized in \cref{fig:channel}(a).
The problem was solved using two overlapping domains as visualized in \cref{fig:channel}(b). The annulus describing the front mesh, visualized in \cref{fig:channel}(b) was chosen to have a width of at least three cells of the background mesh. This is important so that one can identity holes added in the geometry by setting cells to be covered (See \cref{fig:CellTypes}). An interesting effect of this choice is that the front mesh scales with the mesh size, as shown in \cref{Stokesconvmesh}. This means that when we refine the mesh, the number of elements that needs to be deformed does not depend on the total number of degrees of freedom of the MultiMesh.
This means that employing the same deformation scheme on the front mesh would be much more efficient than employing it on the whole domain.

We used the variational formulation of the Stokes equations for two overlapping domains, as derived and discretized in~\cite{johansson2015high},

The shape sensitivity of $J_s$ has been derived in~\cite{schmidtphd} and is
\begin{align}\label{eq:stokesgrad}
  \md J_s(\Omega, u,p)[s]&= \Int{\Gamma_2}{}-s\cdot n
  \left(\der{ u}{n}\cdot\der{ u}{n}\right)\md S.
\end{align}
The shape sensitivity of $J_V$ is obtained by applying the product rule and \cref{HadamardVol}:
\begin{align}\label{eq:volumegrad}
  \md J_V(\Omega)[s]=-2\gamma_1(\vert\Omega\vert-\vert\Omega_0\vert)\Int{\Gamma_2}{}s\cdot n\md S.
\end{align}
Similarly, the shape sensitivity of $J_C$ is obtained using the quotient rule:
\begin{align}\label{eq:centroidgrad}
  \md J_{Cx}(\Omega)[s] &= 2\gamma_2\vert\Omega\vert^{-1}\Int{\Gamma_2}{}s\cdot n(C_x-x)(\text{C}_x-\text{C}_{x0})\md S.
\end{align}
Similar result can be derived for $\md J_{Cy}$.
Combining \cref{eq:stokesgrad,eq:volumegrad,eq:centroidgrad} and obtain the
shape sensitivity
\begin{equation}
  \begin{alignedat}{1}
  \md J(\Omega, u,p)[s]&= \Int{\Gamma_2}{}s\cdot n\Bigg(
   -\left(\der{ u}{n}\cdot\der{ u}{n}\right)
  -2\gamma_1(\vert\Omega\vert-\vert\Omega_0\vert)\\
  &+2\gamma_2\inv{\vert\Omega\vert}\Big[(C_x-x)(\text{C}_x-\text{C}_{x0})  +(C_y-y)(\text{C}_y-\text{C}_{y0})\Big]\Bigg)\md S.
  \end{alignedat}\label{eq:stokes_shape_sensitivity}
\end{equation}
We note that \cref{eq:stokes_shape_sensitivity} does not depend on the adjoint solution.
Hence, for this example one does not have to compute the solution of the adjoint equations.

\subsubsection{Results}\label{sec:stokes_res}
Both a traditional FEM and a MultiMesh FEM Stokes solver was implemented with FEniCS~\cite{alnaes2015fenics,logg2010dolfin}. The meshes used for the experiments are shown in~\cref{Stokesconvmesh}.
Both solvers use a Taylor-Hood element pair for the discretization, i.e. second order piece-wise continuous
polynomials for the velocity and first order piece-wise continuous polynomials for the pressure.
The arising linear systems were solved using the direct solver MUMPS~\cite{MUMPS:1}.
For finer discretizations, it would be beneficial to employ an iterative solver,
though efficient preconditioning for the MultiMesh variational problem has not been properly explored.
The following results use penalty values of $\beta = 6$ in the variational problem defined in~\cite{johansson2015high}, if not otherwise stated.

First, the convergence rates of the Stokes solvers were verified using the manufactured solution given in~\cite{johansson2015high}. For this setup, we expect third and second order convergence rates for velocity and pressure, respectively.
The results, listed in \cref{manufactured_stokes}, show that both solvers achieve the expected convergence-rates.

\begin{figure}[!ht]
  \centering
    \begin{tikzpicture}
      \node [above right, inner sep=0pt]{
        \includegraphics[width=\linewidth]{ManufacturedMeshes.png}
    };
    \node (B) at (0.25\linewidth,-0.25) {(a)};
    \node (B) at (0.75\linewidth,-0.25) {(b)};
  \end{tikzpicture}
    \caption{a) The traditional FEM mesh and b) the MultiMesh FEM meshes for verification of the
      Stokes solver. For the MultiMesh approach, the front mesh is scaled to
        such that the width of the annulus is equivalent to three background cells. In this setup, no extra cells are marked as covered, and there are two interfaces, as shown in \cref{fig:CellTypes}a).}\label{Stokesconvmesh}
\end{figure}
\input{stokes_manufactured}

Next, the shape sensitivity was verified using a Taylor test, as shown in \cref{tab:stokestaylor}. The expected convergence rate, $2$, is obtained. Further tests in random perturbation directions showed similar convergence rates.
We also performed Taylor tests on meshes with different resolution, which revealed that the convergence rate reduces on very coarse meshes due to the discrete inconsistency of the shape sensitivity, see~\cref{sec:shape}.

\input{taylor_stokes.tex}

Finally, we solved the full shape optimization problem.
To ensure that the volume and centroid constraints are sufficiently satisfied,
we solved a sequence of optimization problems, with increasing penalty coefficients $\gamma_1$ and $\gamma_2$.
The optimized mesh of the previous problem is used as an initial mesh for the next optimization problem.
Starting with $\gamma_1=10^5$ and $\gamma_2=10^3$, five optimization problems were solved in which each $\gamma_1$ and $\gamma_2$ were doubled. The number of cut and uncut cells in the initial MultiMesh was 27490 cells.
After a total of 359 optimization iterations, we obtain the final configuration, see~\cref{fig:Stokes90}.
The functional dropped from initially $22.5$ to $18.4$. We observe that the final shape is visually in agreement with~\cite{pironneau1974optimum}, which states that the front and back angle of the object should be $90$ degrees. We measured the front and back angle of our optimized solution to be $87^\circ$ and $84^\circ$ respectively.
The mesh deformation was performed using the advection scheme \cref{eq:ConvEik}. In addition, a centroidal Voronoi tessellation~\cite{du1999centroidal} (CVT) was used on $\Gamma_2$ to preserve the mesh quality near the the front and back wedge.
For the initial mesh, the minimum cell radius ratio was $0.58$, where an equilateral triangle has measure $1$. For the optimized mesh, this had decreased to $0.32$.
Compared to a standard FEM, with a mesh with 25420 cells, an iteration of the optimization algorithm is approximately two thirds with the MultiMesh FEM, compared to standard approach. While the assembly of the linear system for the MultiMesh problem is 4 times slower than for the traditional FEM (1.4 seconds to 0.35 seconds), the deformation of the domain is over 20 times faster (0.35 seconds to 7.85 seconds) for the MultiMesh FEM compared to the traditional FEM. In these timings, we have excluded the execution time of the CVT, as they should be equal for two boundaries with similar resolution. The mesh quality of the optimized traditional FEM solution was 0.37.
\begin{figure}[!ht]
  \centering
  \begin{tikzpicture}
      \node [above right, inner sep=0pt]{
  \includegraphics[width=\linewidth]{Stokes_Optimization.pdf}
     };
    \node (B) at (0.25\linewidth,-0.2) {(a)};
    \node (B) at (0.7\linewidth,-0.2) {(b)};
          \end{tikzpicture}
  \caption{The flow around the initial obstacle. (b) Flow-profile of the fluid around the optimized obstacle. A front and back-angle of 90 degrees are obtained, as proven analytically in Pironneau~\cite{pironneau1974optimum}.}\label{fig:Stokes90}
\end{figure}

\subsection{Orientation of 25 objects in Stokes-flow}
As a final example, we considered the problem of optimally rotate 25 obstacles in Stokes flow to minimize dissipation of energy.
We consider 25 identical objects placed in a structured fashion, as shown in \cref{fig:25stokes}a). There are two identical inlets, with parabolic inlet profiles, and one outlet, that is $66\%$ of total inlet width. The optimization was performed using a MultiMesh consisting of a total of 26 meshes, where each obstacle was represented by a separate mesh, and $22,741$ cut and uncut elements. The stopping criteria of the optimization algorithm was set to  ${(J(\Omega^{k+1})-J(\Omega^k))}/{J(\Omega^{k+1})}<10^{-5}$ and achieved after 50 iterations. The optimized configuration is shown in \cref{fig:25stokes}b).

\begin{figure}[!ht]
  \centering
  \begin{tikzpicture}
      \node [above right, inner sep=0pt]{
        \includegraphics[width=\linewidth]{Stokes_25.pdf}
    };
    \node (B) at (0.25\linewidth,-0.2) {(a)};
    \node (B) at (0.7\linewidth,-0.2) {(b)};
  \end{tikzpicture}
  \footnotesize
  \caption{(a) Initial configuration of 25 objects in a channel with two inlets at the left wall, and one larger outlet at the top. Initial $J_s=72.30$.(b) The optimal orientation of the 25 objects, where each are rotated around their baricenter. Optimal $J_s=56.57$. $50$ iterations where need to obtain the optimized domain. }\label{fig:25stokes}
\end{figure}

config: arxiv