Research Reports of National Institute of Technology, Nagaoka College, Vol. 59, 2024

Paper

Shape optimization of an object in an incompressible
viscous fluid for drag minimization(Considerations for introducing a shape selection
process using deep reinforcement learning)
Yudai SUGIYAMA1, Takahiko KURAHASHI2 and Toshihiko ETO3

Mechanical Engineering, Master’s Program, Graduate School of Engineering, Nagaoka University of Technology

12Department of Mechanical Engineering, Institute of GIGAKU, Nagaoka University of Technology
3

Department of Civil Engineering, National Institute of Technology, Nagaoka College

In this study, we introduce and consider a shape optimization process using deep reinforcement learning
for the drag minimization of an object in an incompressible viscous fluid. Most previous studies on shape
optimization have used the adjoint variable method. However, it has been demonstrated that shape
optimization results are dependent on the initial shape and that the final shape is not necessarily the optimal
shape1). Deep reinforcement learning is a method that learns only from the interactions between an agent
and its environment. We believe that this is appropriate for optimization where the correct answer is
unknown. Furthermore, if the shape design variables increase, the behavior patterns that an agent canconsider will significantly increase. Therefore, we believe that deep learning is more suitable for thispurpose.
Key Words : Deep Q network, Reinforcement Learning, Q-learning, FreeFEM++

1. Introduction

development of machine learning and computer
performance, called deep learning, the field of image
recognition has shown significantly higher
performance than existing methods in general object
recognition contests2). Moreover, a method called
deep Q-learning network (DQN), which learns solely
from trial and error, has been devised and shown tooutperform skilled humans in a majority of games3).
Therefore, this study proposes shape optimization
with a low dependence on the initial shape. A DQN
that learns by searching for solutions to problems for

Shape optimization is a method used in productshape design. Shape analysis with minimum drag andmaximum lift is one of the most important design
processes in shape determination. Shape optimization
using the adjoint variable method has been applied toseveral optimization problems, such as dragminimization problems. However, it has been shownthat shape optimization results depend on the initial
shape and that the final shape is not necessarily the
optimal shape1). However, in recent years, with the6Shape Optimization of an object in an incompressible viscous fluid for drag minimization(Considerations for introducing a shape selection process using deep reinforcement learning)which it does not know the answers is expected towhere, learning rate, action-value function, andperform global optimization.
discount rate are denoted by 𝛼 , 𝑄 and 𝛾 (0 <𝛼, 𝛾 ≤ 1), respectively, and 𝐴 represents the entire
set of possible actions. The second term of Eq.(1)represents the difference between the expected action
2. DQNand the current values.
2.1 What is DQN?
In reinforcement learning (RL), a specific
2.2 Reward
environment is provided, and the behavior is learned
The reward is defined in Eq.(2) as:
such that the reward obtained is maximized. As(2)
𝑟𝑡 = 𝐹𝐷(init) − 𝐹𝐷(𝑡)shown in Fig. 1, RL consists of an agent that makes
behavioral decisions for the model based on thepolicy and an environment that provides rewards.
The environment outputs rewards and the current
state as feedback for the agent’s actions. This agent
learns to perform actions that are expected to be
maximally rewarding in the current state.
RL is commonly used for robot control, Shogi, andGo and is often used for dynamic problems that
depreciate over time. However, other research4),5)
applied it for optimization in static and timeinvariant problems, such as the optimization of thehyperparameters of neural net-works. By receiving
only the value of the first step as feedback in an
inherently time-evolving finite element analysis, theproblem addressed in this study can be considered
static. Various methods have been proposed foraction selection, and this study used a DQN, which
is one such method. The DQN determines actions
based on the action value function, represented by
the following equation. This function was further
modeled using a deep neural network.

Here, the drag forces on the initial shape and the stateat the 𝑡 epoch are denoted by 𝐹𝐷(init) and 𝐹𝐷(𝑡) .
The flow field is an incompressible viscous flow, andthe governing equations are the Navier-Stokes
equations in Eq.(3) and the continuity equation inEq.(4).

𝑢̇𝑖 + 𝑢𝑗 𝑢𝑖,𝑗 + 𝑝,𝑖 −1
𝑢 =0
𝑅𝑒 𝑖,𝑗𝑗

𝑢𝑖,𝑖 = 0

on Ω,

on Ω,

𝑎∈𝐴

𝑢𝑖 = 𝑢̂𝑖 on Γ𝑈 ,𝑇1 = 0, 𝑢2 = 0 on Γ𝑆 ,(1)

Fig. 2 Computational domain and boundary conditions1Fig. 1 Conceptual diagram of RL(4)

where, flow velocity, pressure, and Reynolds
number are denoted by 𝑢𝑖 [m/s], 𝑝[Pa], and 𝑅𝑒,
respectively. The computational domain isconsidered as a typical problem, as shown in Fig.
21). The boundary conditions are as follows:

𝑄(𝑠𝑡 , 𝑎𝑡 ) ← 𝑄(𝑠𝑡 , 𝑎𝑡 )+ 𝛼{𝑟𝑡 + 𝛾 max 𝑄(𝑠𝑡+1 , 𝑎) − 𝑄(𝑠𝑡 , 𝑎𝑡 )}

(3)

4)

7

(5)
(6)

Yudai SUGIYAMA, Takahiko KURAHASHI and Toshihiko ETO
𝑢𝑖 = 0 on Γ𝐵 ,𝑇𝑖 = 0 on Γ𝐷 ,𝑇𝑖 = {−𝑝𝛿𝑖𝑗 + 𝜇𝑢𝑖,𝑗 }𝑛𝑗,

(7)
(8)
(9)

𝑦𝑡𝑄 = 𝑟𝑡 + 𝛾 max 𝑄(𝑠𝑡+1 , 𝑎; 𝜃 − )awhere 𝜃 − are the parameters of the target networkand 𝜃 − = 𝜃𝑡 . Several important changes in network
updates improve the learning stability of the DQN3).
One is that the algorithm is not learned entirely online
but rather by sampling from the experience buffer.
The idea is to use a target network with parameters𝜃 − copied from the main network. The target
network updates its parameters at each defined epoch(𝜃 − ← 𝜃𝑡 ). This method is called a Double DQN(DDQN), which modifies the target value as follows:

where, the stress and unit normal vectors outward
from the boundary are denoted as 𝑡𝑖 , 𝑛𝑗 ,respectively. The fluid force acting on the object isdenoted by 𝐹𝑖 , where the drag and lift forces aredenoted as when 𝑖 = 1, 2 , respectively. 𝐹𝑖 is
obtained using Eq.(10), the integral of 𝑇𝑖 at theboundary Γ𝐵 .
𝐹𝑖 = − ∫ 𝑇𝑖 𝑑Γ.

(14)

(10)

Γ𝐵

Therefore, the reward in Eq.(2) is the difference, that
is, a positive reward if the drag force is smaller than
the initial shape at the 𝑡 epoch. Conversely, if it
becomes larger than the initial shape at the 𝑡 epoch.
it represents a negative reward as a penalty.

𝑦𝑡𝑄 = 𝑟𝑡 + 𝛾 max 𝑄(𝑠𝑡+1 , arg max 𝑄(𝑠𝑡+1 , 𝑎; 𝜃𝑡 ) ; 𝜃 − ) (15)a𝑎where the main network selects the action with thehighest expected action value, as expressed inEq.(15) but adopts the policy 𝜖-greedy method andsometimes randomly selects another action.

2.3 Network update
RL defines the value of an action when 𝑠 and 𝑎are performed according to 𝜋 as follows:

3. Numerical study
𝑄 𝜋 (𝑠, 𝑎) ≔ 𝔼𝑠,𝑎,𝜋 [∑ 𝑟𝑡 ](11)

Shape optimization using the DQN and adjoint
variable method was performed based on a previous
study’s analytical conditions1). The flow field wasanalyzed using FreeFEM++ v4.9, and the free
software application focused on solving partial
differential equations using the finite element method.
The programming language Python 3.8.0 was used to
implement the DQN, and the machine learning
framework PyTorch 1.13.1 was used.

where the policy for choosing action is denoted by 𝜋.
The optimal values are as follows:𝑄 ∗ (𝑠, 𝑎) ≔ max 𝑄 𝜋 (𝑠, 𝑎)𝜋(12)Instead of using tabular encodings, neural networks
are used to learn estimates by parameterized𝑄(𝑠, 𝑎; 𝜃) to scale up to large problems such asmodels with many design variables. The following
formula is used to update the network parameters.

3.1 Analysis conditions
Fig. 3 and 4 show the finite element mesh, andcomputational domain, respectively. The parameters
used in the DQN analysis are listed in Tab. 1. In the𝜃𝑡+1 ← 𝜃𝑡 + 𝛼 (𝑦𝑡𝑄 − 𝑄(𝑠𝑡 , 𝑎𝑡 ; 𝜃𝑡 )) ∇𝜃 𝑄(𝑠𝑡 , 𝑎𝑡 ; 𝜃𝑡 ) (13)finite element analysis, the time step width was
0.03[s], end time step was 2,000[s], and the Reynolds
𝑄
number was 100.
where、𝑦𝑡 is target value as follows:8

Shape Optimization of an object in an incompressible viscous fluid for drag minimization(Considerations for introducing a shape selection process using deep reinforcement learning)where, the design variable in the DQN analysis is𝑎 and the area of the initial shape was determined
given by the elliptical equation as in Eq.(16).
using Newton’s method. The Adam optimization
method, which is a sequential optimization method22that uses first-order derivatives of the objective
𝑥
𝑦
(16)
+ 2=1
2
function with various mod-ifications to the steepest𝑎𝑏descent method, was used. Specifically, this method
combines RMSprop and mom-entum SGD. Eq.(1) is
The action was to increase or decrease 𝑎. To update
represented by a neural network, in which the second
the shape with a constant area, 𝑏 was obtained fromterm is the difference between the expected action
value and the current expected value, as follows:
1

𝑁

𝐸 = 2∑𝑄 2(𝑄 𝑛 (𝑠𝑡 , 𝑎; 𝜃) − 𝑦𝑡 )

(17)

𝜕𝐸

(18)

𝑛=1

∇𝐸 ≔ 𝜕𝑤

where, 𝐸 and ∇𝐸 in Eqs.(17) and (18) denote theerror function and gradient of the error function,
respectively. By applying Adam to them we obtain:

Fig. 3 Finite element mesh1)

[𝑡]

[𝑡−1]

[𝑡]

[𝑡−1]

𝑚𝑖 = 𝛽1 𝑚𝑖𝑣𝑖 = 𝛽2 𝑣𝑖

Fig. 4 Computational domain1)

Max epoch
Max EpisodeWorm up epoch
Batch size
Learning rate
Adam parameter 1
Adam parameter 2
Adam parameter 3
EpsilonGamma
Initial design variable
Action width
Newton threshold

𝛼
𝛽1
𝛽2
𝜖1
𝜖2
𝛾
a, 𝑏
𝜖3

[𝑡]

200
100
10
5
0.00025
0.9
0.999

[𝑡]

where, 𝑚𝑖+ (1 − 𝛽2 )𝛻(𝐸 [𝑡] )2

(20)

𝑚𝑖

(21)

(1−𝛽1𝑡 )
[𝑡]

𝑣𝑖

[𝑡]

𝑣̂𝑖 =
𝜃𝑖

(19)

[𝑡]

[𝑡]

𝑚
̂𝑖 =

Tab.1 DQN analysis conditions+ (1 − 𝛽1 )𝛻𝐸(22)(1−𝛽2𝑡 )
[𝑡−1]

= 𝜃𝑖

[𝑡]

−𝛼

̂𝑖
𝑚

(23)

[𝑡]

(√𝑣̂𝑖 +𝜖1 )in Eq.(19) is the momentum forsuppressing the ith parameter vibration at 𝑡, 𝑚0 =0 . 𝛽1 is the decay rate, the hyperparameter that
determines how much of the past gradient is

1× 10−80.1
0.9
0.5[m]
0.01[m]accumulated, and

[𝑡]

𝑣𝑖

in Eq.(20) denotes theaccumulation of past gradients for the ith learningparameter at time 𝑡. The initial value was zero. The
ratio of the past gradient to the latest gradientaccumulated was 𝛽2 ∶ 1 − 𝛽2 . Therefore, past

1 × 10−6

9Yudai SUGIYAMA, Takahiko KURAHASHI and Toshihiko ETO
information accumulates exponentially, allowing
more past information to be largely ignored as theepochs progressed, and eliminating the gradient loss
problem.

Fig. 5 Transition from initial to convergent

smooth shape.
Immediately after the inflow, the drag was large,
with drags of 39.88 and 29.67 for the DQN andadjoint variable methods at the first step, respectively.
Fig. 8 shows the drag after the second step. The
cumulative values of the square of the drag for Steps
1 to 2,000 were 26.42 and 50.46 for the DQN andadjoint variable methods, respectively. It is difficult
to represent the complexity of the geometry in thisanalysis because of the small number of design

3.2 Result
Fig. 5 shows the transition from the initial toconvergent solutions during the verification of theinitial shape, where 𝑎 is the width from the center ofthe ellipse, indicating gradual widening andconvergence to a width of 2.14 for an initial conditionof width 0.5. The width of the movement in one
action was 0.01[m]. Fig. 6 shows the final shapes ofshape optimization using DQN and the associatedvariable method.
The adjoint variable method seeks a bullet shapewith a pointed tip and the rounded rear. The pressure
distributions after 2,000 steps for the two geometries
are shown in Fig. 7. The drag force history for 2,000
steps is shown in Fig. 8. The final shape of theadjoining variable method appears to be a backward
vortex, and the flow appears to peel off around the
middle of the object. Pressure resistance is dominant
in vortex-generating objects. Contrarily, there were
no noticeable eddies in the DQN shape, suggestingthat almost no exfoliation occurred. Therefore, we
believe that the frictional drag was dominant, but thefrictional drag was not that large because of theFig. 6 Final shape of DQN (top) andadjoint variable method (bottom)

800

Reward [-]600

𝑎 = 2.14
400

𝑎 = 0.5
200

0
0.1

1

10

𝑎 [m]

10

Shape Optimization of an object in an incompressible viscous fluid for drag minimization(Considerations for introducing a shape selection process using deep reinforcement learning)
Fig. 7 Pressure distribution for DQN (top) andAcknowledgement: We would like to thankadjoint variable method (bottom)
Editage (www.editage.com) for English language
editing.
variables in the DQN. However, we believe that
resetting the shape and repeating the search for each
episode is less likely to result in a locally optimal
References1)
Hiroko Yagi and Mutsuto Kawahara: Optimal shapesolution. We also believe that a larger step width anddetermination of a body located in incompressible viscouswider search compared with conventional analysis
fluid flow, Computer Methods in Applied Mechanics andconditions contributed to the success of the search.
Engineering, 196, 5084-5091, 1 November 2007.
0.4

2)

Simonyan, K. and Zisserman, A:Very Deep Convolutional
Networks for Large-Scale Image Recognition. arXiv

0.35

preprint arXiv:1409.1556, 2014.
3)

Drag Force [N]0.3

M.G. Bellemare, A. Graves, M. Riedmiller, A.K. Fidjeland,

DQN
Adjoint

0.25

V. Mnih, K. Kavukcuoglu, D. Silver, A.A. Rusu, J. Veness,

G. Ostrovski, S. Petersen, C. Beattie, A. Sadik, I.
Antonoglou, H. King, D. Kumaran, D. Wierstra, S. Leggand D. Hassabis: Human-level control through deep

0.2

reinforcement learning, Nature., 518, 529–533, 2015.
0.15

4)

S. S. Richard and G. B. Andrew: Reinforcement Learning:
An Introduction, MIT Press., 2014, 2015.

0.1

5)

K. Li, J. Malik: Learning to optimize. ICLR 2017conference., 5 November 2016 (modified: 22 July 2022).

0.05
6)

M. Andrychowicz, M. Denil, S. Gomez, M. W. Hoffman,
D. Pfau, T. Schaul, B. Shillingford and Nando de Freitas:

0

1

10

100

1000

Learning to learn by gradient descent by gradient descent,Time step [-]arXiv preprint arXiv, 1606, 044-74v2, 14 Jun 2016.

Fig. 8 Drag history(Received July 19, 2023)

4. Conclusion
In this study, we proposed a shape optimization
method using a DQN and compared it with aconventional shape optimization method that uses
adjoint variable method. Although the number ofdesign variables was different and therefore, notcompletely comparable, we believe that the DQN is
useful for shape optimization of models with fe w
design variables based on drag forces. In the future,
we intend to verify this using additional design
variables.

11