# Provably Safe Deep Reinforcement Learning for Robotic Manipulation in Human
Environments

Jakob Thumm and Matthias Althoff The authors are with the Department of
Informatics, Technical University of Munich, 85748 Garching, Germany.
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Deep reinforcement learning (RL) has shown promising results in the motion
planning of manipulators. However, no method guarantees the safety of highly
dynamic obstacles, such as humans, in RL-based manipulator control. This lack
of formal safety assurances prevents the application of RL for manipulators in
real-world human environments. Therefore, we propose a shielding mechanism
that ensures ISO-verified human safety while training and deploying RL
algorithms on manipulators. We utilize a fast reachability analysis of humans
and manipulators to guarantee that the manipulator comes to a complete stop
before a human is within its range. Our proposed method guarantees safety and
significantly improves the RL performance by preventing episode-ending
collisions. We demonstrate the performance of our proposed method in
simulation using human motion capture data.

##  I INTRODUCTION

In recent years, researchers solved many complex manipulation tasks using deep
reinforcement learning (RL), such as operating door handles [1], playing table
tennis [2], stacking boxes [3], and controlling multiple robotic arms [4].
Furthermore, recent work [5, 6, 7] demonstrated that RL-controlled
manipulators could successfully maneuver in environments with dynamic
obstacles, dodging moving obstacles, and reaching the goal state consistently.
Despite these promising results, one open challenge in RL is to formally
guarantee the safety of surrounding humans due to their unpredictable movement
and many degrees of freedom (DoF). The lack of safety assurances makes
applying any vanilla RL method in human working environments irresponsible.

To overcome these restrictions, we propose a safety shield concept for RL that
guarantees human safety at all times, as shown in Fig. 1. In essence, our
safety concept is speed and separation monitoring (SSM) according to DIN EN
ISO 10218–1 2021, 5.10.3 [8] warranting that the robot is in a completely
stopped state before any collision with a human could occur. We need to verify
the robot trajectory at a high frequency to achieve this strong safety
criterion while still being able to maneuver close to humans. However, for the
RL agent it suffices to output actions at a low frequency to make long-term
decisions like moving around a human to reach its goal. Our safety shield
combines the low-frequency RL agent with our high-frequency formal
verification. The RL agent used in this work is a state-of-the-art soft actor-
critic (SAC) [9] with hindsight experience replay (HER) [10]. However, our
safety shield can be used with any online RL agent. This work presents the
first provably safe robot manipulator control based on deep RL in human
environments.

Figure 1: Structure of the proposed RL agent with safety shield.

###  I-A Related work

Gu et al. [1] were one of the first to show that various complex manipulation
tasks can be learned using only off-policy deep RL methods. They presented a
method to train multiple real robots in parallel to solve a challenging door-
opening task. Shortly later, [10] proposed HER to achieve randomized
manipulation goals. The combination of HER with off-policy RL contributed to
many remarkable results in pushing, pick and place, and throwing tasks as
presented in [10, 11, 12]. Recently, HER was combined with state-of-the-art RL
algorithms like SAC to solve complex manipulation tasks such as controlling
multiarm manipulators [4].

There have been many approaches that introduce safety constraints to the
exploration process of RL agents. Simple reward shaping has shown to be
insufficient to guarantee safety as it is infeasible for universally deciding
between fulfilling the safety constraints and rapidly accomplishing the goal
[13, 14]. Many recent works use constrained policy optimization techniques
[13, 15, 14] to reduce the number of safety-critical interactions with the
environment. Although some of these approaches show a very low number of
safety constraint violations, none can warrant safety at all times.

In many scenarios, especially when interacting with humans, a very low chance
of collision is still unacceptable. Therefore, recent works [6, 16, 17, 18] in
the field of provably safe RL (often referred to as shielding, safety layers,
sandboxing, or constraining bolts) focus on provable safety guarantees. The
core idea of these methods is to ensure that only safe actions are sent to the
environment. As a result, unsafe actions are either replaced with safe
actions, projected to the safe action space, or entirely prevented by limiting
the action space. Many provably safe RL approaches have been proposed for
simple environments such as Atari games [19], grid worlds [20, 21], or board
games [22]. Most notably, Hunt et al. [17] exhibited how safe end-to-end
learning from image data could be achieved in discrete action spaces.
Unfortunately, these methods rely on a deterministic environment to assure
safety, so they are not directly applicable to a complex task like human–robot
collaboration or coexistence. For the use case of RL for manipulators, Pham et
al. [6] proposed a differentiable safety layer (OptLayer) that projects any
unsafe action to the closest safe action that satisfies the given constraints.
Despite promising results, this method comes with two flaws. First, the
OptLayer is only formally correct if the constraints of the underlying
quadratic program are precisely fulfilled, and the optimizer finds an exact
solution in the given time. Second, it is difficult to formalize the safety
constraints for unknown and complex human motions. The reachable set-based
safety layer of [16] provides formal safety for autonomous driving highway
scenarios by masking out all nonsafe actions. Training an RL agent with this
safety layer leads to faster convergence while all collisions can be avoided.
This method was designed for a discrete action space and is therefore not
directly transferable to the manipulator application. The approach closest to
ours is the reachability-based trajectory safeguard of [18]. They use
reachability analysis to predict if the agent could collide with static
obstacles on the current trajectory and replan a safe trajectory if necessary.
However, this approach is only designed and tested for static environments.

###  I-B Contributions

We present a novel safety shield111Our code and models are publicly available
at https://github.com/JakobThumm/safe_rl_manipulators. that replaces unsafe
actions from the agent with provably safe actions from a high-frequency safe
trajectory planner (failsafe planner). To the best of our knowledge, our
method is the first safety measure for RL that provides provable safety for
continuous action spaces in high-dimensional state-spaces and hard-to-predict
dynamic environments. Our safety shield can easily be applied to a large
variety of manipulation tasks. Compared to all previously mentioned provably
safe RL approaches, we sample from the entire set of actions (instead of a
safe subset) and check the safety of the current action during execution at
high frequency. This allows for quick reactions to highly dynamic human
motion. To summarize, our safety shield guarantees that the robot stops before
a collision with a human could occur while still allowing the maximal freedom
of movement under the given safety constraints.

###  I-C Article structure

Section II summarizes our problem and presents the RL basics and notation used
in the article. Next, Section III presents our proposed safety shield and RL
agent. We then discuss the experimental setup and results in Section IV and
present our conclusion in Section V.

##  II PRELIMINARIES

###  II-A Problem statement

In our RL setting, a six DoF modular robot, as described in [23], is mounted
on a working table, where a human repeatedly performs a task. In each episode,
the agent has to reach a randomized goal joint position, further referred to
as episode goal $g$, from a fixed initial position and evade the human. We
only consider episode goals that are collision-free with the static
environment. This paper aims to provide formal safety guarantees for nearby
humans with unknown motion behavior while consistently reaching the goal.
Hereby, the RL agent receives the following observations: the current joint
position and velocity, episode goal, current Cartesian end-effector position,
and the relative Cartesian positions of the human wrists and head in relation
to the end-effector. Since the only static obstacle in our scene is the table,
we do not add static obstacle information to the observation space. Fig. 2
shows our setup and illustrates the observations.

Figure 2: Simulation setup with the robot mounted on a desk and the human
performing a task. The provided observations are listed on the right hand
side. Figure 3: Illustration of intended and failsafe trajectories. In time
steps $t_{0}\dots t_{2}$, the intermediate goal for the intended trajectories
is $\bar{g}_{i}$. Before time step $t_{0}$, a failsafe trajectory is
calculated for time step $t_{1}$. In between $t_{0}$ and $t_{1}$, an intended
trajectory for $t_{1}$ followed by a failsafe trajectory for $t_{2}$ is
calculated. The intended trajectory in $t_{1}$ is not safe, therefore the
robot executes the failsafe trajectory. A new intended trajectory for $t_{2}$
followed by a failsafe trajectory for $t_{3}$ is calculated during execution
of the failsafe trajectory between $t_{1}$ and $t_{2}$. At $t_{2}$, the
intended trajectory is verified as safe and can be followed. Before the
trajectory calculation for $t_{3}$, the agent outputs a new intermediate goal
$\bar{g}_{i+1}$, so that the new intended trajectory is calculated from
$t_{3}$ to $\bar{g}_{i+1}$. Since the intended trajectory to $\bar{g}_{i+1}$
is not verified as safe at $t_{3}$, the robot follows the previous failsafe
trajectory. A new intended trajectory from $t_{4}$ to $\bar{g}_{i+1}$ is
calulated followed by a failsafe trajectory. This intended trajectory is safe,
so that the robot switches to it.

###  II-B RL

For RL, we use the notation in [9] and consider the Markov decision process
defined by the tuple $\left(\mathcal{S},\mathcal{A},p,r\right)$ with both
continuous state space $\mathcal{S}$ and action space $\mathcal{A}$. For
simplicity, we assume the state to be fully observable. The transition
function
$p:\mathcal{S}\times\mathcal{S}\times\mathcal{A}\rightarrow\mathbb{R}$ denotes
the probability density function of reaching the next state
$\boldsymbol{s}_{i+1}$ when choosing action $\boldsymbol{a}_{i}$ in state
$\boldsymbol{s}_{i}$. After each transition, the agent receives a reward from
the environment according to the reward function
$r:\mathcal{S}\times\mathcal{A}\rightarrow\mathbb{R}$. The agent learns a
stochastic policy $\pi(\boldsymbol{a}_{i}|\boldsymbol{s}_{i})$ for action
$\boldsymbol{a}_{i}$ given state $\boldsymbol{s}_{i}$.

The optimal policy is learned using the online off-policy SAC algorithm first
presented in [9]. Our specific implementation extends the spinningup [24]
version of SAC with HER and is further described in Section III-B. The goal of
SAC is to automatically balance exploration and exploitation by maximizing the
accumulated reward and the information content of the policy function.
Therefore, the entropy $\mathcal{H}$ (see e.g. [25, Eq. 4.2]) of the policy is
used as part of the optimization objective (from [9, Eq. 1])

$J(\pi)=\sum_{t=0}^{T}\mathbb{E}_{\left(\boldsymbol{s}_{i},\boldsymbol{a}_{i}\right)\sim\rho_{\pi}}\left[r\left(\boldsymbol{s}_{i},\boldsymbol{a}_{i}\right)+\alpha\mathcal{H}\left(\pi\left(\cdot\mid\boldsymbol{s}_{i}\right)\right)\right].$
(1)

In many robotic applications, one encounters the problem of sparse rewards,
where the agent needs to get to a random and hard-to-reach goal. If the agent
rarely reaches the goal, almost no transitions with goal encounters are added
to the replay buffer. In order to create more of these goal transitions, [10]
proposed HER with the idea to sample additional transitions with fictional
goals, where the goal $g$ of an episode is added to state
$\boldsymbol{s}_{i}$, written as $\boldsymbol{s}_{i}||g$. The goal-conditioned
transition tuples can be written as
$\left(\boldsymbol{s}_{i}||g,\boldsymbol{a}_{i},r_{i},\boldsymbol{s}_{i+1}||g\right)$.
For each transition, additional fictional episode goals
$\tilde{g}^{h},\,h=1\dots k_{\text{HER}}$ are created from states that were
actually reached in the episode of the transition, leading to $k_{\text{HER}}$
(in our case $k_{\text{HER}}=4$) new transitions
$\left(\boldsymbol{s}_{i}||\tilde{g}^{h},\boldsymbol{a}_{i},\tilde{r}_{i}^{h},\boldsymbol{s}_{i+1}||\tilde{g}^{h}\right)$.
By sampling new goals from states that have actually been reached, many
successful goal transitions are added. This method enables us to use a sparse
reward, where $r_{i}=0$ if the goal was reached and $r_{i}=-1$ otherwise.
Since we terminate an episode if a collision occurs, e.g., the human walks
into the stationary robot, a negative reward for collisions is not strictly
necessary.

##  III METHODOLOGY

###  III-A Safety shield

Our proposed safety shield is based on the provably safe trajectory planner,
called failsafe planner in [26]. The shield functionality is illustrated in
Fig. 1 and further elaborated in this subsection. Please refer to [26, 27, 23,
28] for a more detailed description of the trajectory planning and
verification process.

In our application, the actions of the RL agent are intermediate goal joint
positions $\bar{g}_{i}$. These intermediate goals $\bar{g}_{i}$ are not to be
confused with the episode goal $g$ or the fictional episode goals for HER
$\tilde{g}^{h}$. In each episode, the agent outputs multiple intermediate
goals $\bar{g}_{i}$ to reach the single episode goal $g$. The HER goals
$\tilde{g}^{h}$ are only relevant for training the agent. Furthermore, we
would like to emphasize that the RL agent outputs the intermediate goals at a
lower frequency than our safety shield is operated. Each RL action is executed
for a user-defined time interval $\Delta T$, or until all joints are within a
small $\epsilon$-range of the intermediate goal. We refer to this as one RL
step and use the index $i$. The user-defined time between two safety shield
verification steps is $\Delta t$ with $\Delta t<\Delta T$ and the index $k$ is
used.

#### Long-term and failsafe planning

This paragraph describes the trajectory planning in our safety shield and is
accompanied by a simplified example in Fig. 3. To reach the intermediate goal
$\bar{g}_{i}$, we calculate an intended trajectory from the current joint
state $\boldsymbol{x}_{k}$ to $\bar{g}_{i}$ using a long-term planner. We
guarantee that a failsafe trajectory always exists, which brings the robot to
an inherently safe state, i.e., a full stop, because a dynamic obstacle can
block the robot’s path during the execution of this intended trajectory.
Formal safety can now be guaranteed by induction, where it is assumed that the
robot starts in a safe state (stopped). The first intended trajectory is
executed starting from $t_{0}$ if no collision with a human is possible
between $t_{0}$ and $t_{1}$ and if a collision-free failsafe trajectory exists
starting from $t_{1}$. During the execution of a trajectory between $t_{k}$
and $t_{k+1}$, we calculate a failsafe trajectory starting from $t_{k+2}$. If
this failsafe trajectory cannot be verified as safe, we execute the previously
verified failsafe trajectory starting from $t_{k+1}$. Since each failsafe
trajectory ends in a safe state, we can guarantee safety for an infinite time
horizon. Similarly, one time step of the intended trajectory between $t_{k}$
and $t_{k+1}$ is only executed if it is verified as safe; otherwise, the
failsafe trajectory is executed. In each time step, a new intended and
failsafe trajectory is planned.

In our RL scenario, the agent repeatedly updates the intermediate goal. If
there is a new intermediate goal $\bar{g}_{i+1}$ at time step $t_{k}$, the
intended trajectory is calculated from the state at time $t_{k+1}$ to
$\bar{g}_{i+1}$, and a new failsafe trajectory is calculated from the new
intended trajectory starting at time $t_{k+2}$. If the new intended and
failsafe trajectories can be verified as safe, the robot changes to the new
intended trajectory. Otherwise, it follows the last verified failsafe
trajectory at $t_{k+1}$.

#### Adaption to manipulator control

In our use case, the 3D trajectory planning can become computationally
expensive, especially in regards to the high-frequency requirement. Therefore,
we do not recompute the intended trajectory in every time step. Instead, we
only compute intended trajectories for every new intermediate goal
$\bar{g}_{i+1}$ until an executable one is found and use velocity scaling for
the failsafe trajectory as described in [26]. In short, a failsafe trajectory
is path-consistent with the intended trajectory and ends in the complete stop
of the robot. We use the synchronous Type IV online trajectory planner of [29]
for the intended and the failsafe planner. Thus, the intended trajectory
planning adheres to velocity, acceleration, and jerk limits
($v^{\text{traj}}_{\text{max}}$, $a^{\text{traj}}_{\text{max}}$, and
$j^{\text{traj}}_{\text{max}}$), which are gentle for the manipulator joints.
To achieve fast-braking, we increase the acceleration and jerk limits for the
failsafe planning ($a^{\text{failsafe}}_{\text{max}}$ and
$j^{\text{failsafe}}_{\text{max}}$ ) to the physical limits given by the robot
manufacturer.

#### Verification

The trajectory verification is based on collision checking of the reachable
sets of all possible human motions and the intended robot trajectory. It is
assumed that no human joint moves faster than
$2\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$ as defined in DIN EN ISO 13855:2010
[30], and all human joint positions are reliably measurable within a specified
error bound, e.g., by using a motion capture system. The safety verification
also works with other perception methods like light curtains, as presented in
[31], but the robot’s movement would become more conservative. The entire
space a human can occupy within the time interval $\left[t_{k},t_{k+b}\right]$
is defined as the reachable occupancy
$\Gamma\left(\left[t_{k},t_{k+b}\right]\right)$, where $b$ is the number of
time steps needed to bring the robot to a complete stop, and is calculated
using the task space approach presented in [27, 23].

The human and robot occupancies are modeled with capsules, as shown in Fig. 4,
to achieve fast computation times. A capsule comprises a cylinder with half-
spheres at both ends and is defined by a line segment $l({p}_{1},{p}_{2})$
with endpoints ${p}_{1}$ and ${p}_{2}$ and a radius $r$. The robot occupancy
is described as a set of capsules enclosing each robot link’s movement between
time $t_{k}$ and $t_{k+b}$. Hereby, we follow the approach described in [26]
to guarantee over approximitaion of the link occupancies. This approach
assumes that the controller follows the desired trajectory exactly. If
additionally deviation occurring due to low-level control should be
considered, the approach presented in [31] can be used. A robot motion is
verified as safe if no robot capsule intersects with any human capsule for all
times before reaching the resting position.

Figure 4: Example of a human–robot capsule pair in collision.

###  III-B Safe RL

1 Initialize:

2

  * •

soft actor-critic agent $\mathbf{A}$

  * •

environment $\mathbf{E}$

  * •

replay buffer $R$ and local buffer $L$

  * •

$T_{\text{total}}\leftarrow 0$, $T_{\text{last\\_update}}\leftarrow 0$

  * •

randomly initialize weights of $\pi_{b}$

3for $j\leftarrow 0$ to $n_{\text{epochs}}$ do

4 for $l\leftarrow 0$ to $n_{\text{episodes\\_per\\_epoch}}$ do

5 reset $\mathbf{E}$ and clear $L$ ;

6 for $i\leftarrow 0$ to $t_{\text{max\\_episode}}$ do

7 if $T_{\text{total}}\geq k_{\text{start\\_steps}}$ then

8 $a_{i}\leftarrow\pi_{b}\left(s_{i}||g\right)$;

9

10 else

11 $a_{i}\leftarrow$ random action;

12

13 $\left(s_{i+1},\,r_{i},\,done\right)\leftarrow\mathbf{E}$.step($a_{i}$) ;

14
$L\left[i\right]\leftarrow\left(s_{i}||g,\,a_{i},\,r_{i},\,s_{i+1}||g\right)$;

15 $T_{\text{total}}\leftarrow T_{\text{total}}+1$;

16 if done then break ;

17

18 $l_{episode}\leftarrow i$;

19 for $i\leftarrow 0$ to $l_{episode}-1$ do

20 $\left(s_{i}||g,\,a_{i},\,r_{i},\,s_{i+1}||g\right)\leftarrow
L\left[i\right]$;

21 $R$.store$\left(s_{i}||g,\,a_{i},\,r_{i},\,s_{i+1}||g\right)$;

22 remaining_ids $\leftarrow(i+1)\dots(l_{episode}-1)$;

23 randomly select $k_{\text{HER}}$ random_ids from remaining_ids;

24 for id in random_ids do

25 $\left(s_{id}||g,\,\dots\right)\leftarrow L\left[id\right]$;

26 $\tilde{g}\leftarrow s_{id}$;

27 $\tilde{r}_{i}\leftarrow\mathbf{E}$.ComputeReward($s_{i+1}||\tilde{g}$);

28
$R$.store$\left(s_{i}||\tilde{g},\,a_{i},\,\tilde{r}_{i},\,s_{i+1}||\tilde{g}\right)$
;

29

30 if $T_{\text{total}}\geq T_{\text{update\\_after}}$ and
$(T_{\text{total}}-T_{\text{last\\_update}})\geq T_{\text{update\\_every}}$
then

31 for $m\leftarrow 0$ to $T_{\text{total}}-T_{\text{last\\_update}}$ do

32 $\mathbf{A}$.Update($R$.sample_batch());

33

34 $T_{\text{last\\_update}}\leftarrow T_{\text{total}}$;

35

Algorithm 1 Training procedure

Our policy network has a $\tanh$ activation function for the output layer and
thus outputs actions
$a_{i}=[a^{0}_{i},a^{1}_{i},\dots,a^{N}_{i}],a^{n}_{i}\in[-1;1]\forall n\in
1\dots N$ with $N$ being the number of joints. We convert these actions to
intermediate goal joint positions with
$q^{n}_{\bar{g},i+1}=q^{n}_{i}+a^{n}_{i}\Delta q_{\text{max}}$, where
$q^{n}_{i}$ is the current position of joint $n$ and $\Delta q_{\text{max}}$
is the maximum joint position difference per RL step. In recent literature, it
is more common to use the relative change in 3D position of the robot end-
effector as action space and calculate $q^{n}_{i+1}$ using inverse kinematics
(e.g. in [3]) since it leads to a lower dimensional action space. However, we
use joint positions as actions to have full control over each joint.
Nevertheless, both approaches are equally compatible with our proposed safety
shield. The choice of $\Delta q_{\text{max}}$ depends on the maximum velocity
of the robot $v^{\text{traj}}_{\text{max}}$ and the execution time of each RL
step $\Delta T$. Choosing $\Delta q_{\text{max}}\geq
v^{\text{traj}}_{\text{max}}\Delta T$ allows the agent to output a new action
before the current action is finished. This is desirable because each
trajectory ends in a stopped state and we do not want the robot to stop after
each action. If the robot collides with the static environment in its
intermediate goal, a new action is generated randomly with a uniform
distribution over the entire action space until a collision-free intermediate
goal is determined. An alternative to this method is to project the unsafe
action to the safe action space with respect to the static environment.
However, this is nontrivial in the case of a high-dimensional manipulator.

In our training procedure, which is described in Algorithm 1, an episode is
done if the episode goal is reached, the robot is in collision, or by timeout
when the maximum number of RL steps $T_{\text{max\\_episode}}$ is reached. We
consider the episode goal as reached, if
$|q^{n}_{i+1}-q^{n}_{g}|<\epsilon_{g},\forall n\in 1\dots N$, where
$q^{n}_{i+1}$ is the position of joint $n$ after the environment step,
$q^{n}_{g}$ is the episode goal position of joint $n$, and $\epsilon_{g}$ is a
small user-defined tolerance for reaching the goal. At the beginning of each
episode, the environment is reset (Algorithm 1 Algorithm 1) by moving the
robot back to its initial position and resetting the human animation. The
action of a single RL step is executed in Algorithm 1 by the function
$\mathbf{E}$.step($a_{i}$). To achieve a high exploration rate in the
beginning of the training, we choose random actions instead of the policy for
the first $k_{\text{start\\_steps}}$ RL steps (Algorithm 1). After each
episode, HER samples are added to the replay buffer as described in Algorithm
1 to Algorithm 1. The agent networks are first updated after
$T_{\text{update\\_after}}$ RL steps to improve exploration at the start of
training and then every $T_{\text{update\\_every}}$ RL steps from Algorithm 1
to Algorithm 1. All relevant hyperparameters for RL training and operation are
given in Table I.

TABLE I: Hyperparameters

Parameter | Value
---|---
$\Delta t$ | $4\text{\,}\mathrm{m}\mathrm{s}$
$\Delta T$ | $200\text{\,}\mathrm{m}\mathrm{s}$
$v^{\text{traj}}_{\text{max}}$ | $2\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}\mathrm{/}\mathrm{s}$
$a^{\text{traj}}_{\text{max}}$ | $2\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}\mathrm{/}\mathrm{s}^{2}$
$j^{\text{traj}}_{\text{max}}$ | $15\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}\mathrm{/}\mathrm{s}^{3}$
$a^{\text{failsafe}}_{\text{max}}$ | $10\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}\mathrm{/}\mathrm{s}^{2}$
$j^{\text{failsafe}}_{\text{max}}$ | $400\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}\mathrm{/}\mathrm{s}^{3}$
$\Delta q_{\text{max}}$ | $0.4\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}$

Parameter | Value
---|---
Hidden layers (all networks) | $3$
Hidden units per layer | $64$
Minibatch size | $128$
Replay buffer size | $10^{6}$
Discount factor $\gamma$ | $0.99$
Entropy tradeoff $\alpha$ | $0.2$
$k_{\text{start\\_steps}}$ | $5000$
$T_{\text{update\\_after}}$ | $1000$
$T_{\text{update\\_every}}$ | $200$
$n_{\text{epochs}}$ | 200
$n_{\text{episodes\\_per\\_epoch}}$ | 30

[trim=0.0cm 0.0cm 0.0cm 0.0cm, clip, mode=buildnew,
width=]tikz/randomized_goal_hor

[trim=0.0cm 0.0cm 0.0cm 0.0cm, clip, mode=buildnew,
width=]tikz/evade_humans_hor

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/collision_rate_6dof_rand

(a)

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/done_rate_6dof_rand

(b)

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/episode_ending_safe

(c)

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/collision_rate_6dof_fixed

(d)

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/done_rate_6dof_fixed

(e)

[trim=0.2cm 0.2cm 0.2cm 0.0cm, clip, mode=buildnew,
width=0.3]tikz/episode_ending_safe_fixed

(f)

Figure 5: Evaluation of the randomized-goal (upper row) and human-evasion
(lower row) experiment. Figures (a) and (d) compare the rate of safety-
critical collisions per episode of the baseline and the safe agent. The
success rate of both experiments is displayed in (b) and (e). The reason for
episode endings using the safe agent are displayed in (c) and (f). Hereby, a
safe collision occurs when the human touches the stopped robot.

###  III-C Software structure

Our software framework is based on Gazebo and ROS Noetic. The robot simulation
uses the Open Dynamics Engine (ODE) solver with a step size of
$1\text{\,}\mathrm{m}\mathrm{s}$. We choose the ODE QuickStep method with $10$
solver iterations. The safety shield is implemented in C ++ 222Part of our
code is based on the open source reachability analysis tool SaRA (available at
https://github.com/Sven-Schepp/SaRA) [32]. to achieve fast calculation times
with an average of $0.5\text{\,}\mathrm{m}\mathrm{s}$ per verification
step333Run on a ThinkPad P15 Gen 1 with an Intel(R) Core(TM) i7-10750H CPU and
$32\text{\,}\mathrm{G}\mathrm{B}$ DDR4 RAM.. With the safety shield operating
at $250\text{\,}\mathrm{H}\mathrm{z}$, our simulation runs roughly five times
faster than real-time.

##  IV RESULTS

###  IV-A Experimental setup

We conduct two different experiments. In the first randomized-goal experiment,
the episode goal is randomly and uniformly sampled over the entire joint
space, and the start state is fixed with all joints at their zero position.
Notably, the human may block the episode goal. For the second human-evasion
experiment, we manually define the start and episode goal position so that the
robot collides with the human if it takes the shortest path. The agent must
learn to evade the human by taking a longer path or waiting until the human
moves away from the table. In this scenario, the episode goal is only slightly
randomized so that the human never blocks the goal state. The two experiments
are also shown in our accompanying video444Available at
https://youtu.be/Lzrs2HQUIOc.

For the simulation, the human motion is taken from real-world motion capture
data provided by CMU555Data publicly available at http://mocap.cs.cmu.edu/.
Currently, we are using animation 62_01, containing a human walking to the
working table, where they perform a wrenching action. This motion pattern is
diverse and complex and can effectively validate our concept. To prevent
overfitting to the motion data, we randomize the $\mathsf{x}$\- and
$\mathsf{y}$-position and start time of the human animation in each episode
uniformly in the range of $-0.2\text{\,}\mathrm{m}0.2\text{\,}\mathrm{m}$ and
$0\text{\,}\mathrm{s}1\text{\,}\mathrm{s}$.

In all experiments, the maximum episode length is set to $100$ RL steps per
episode to give the robot adequate time to move around the human to any goal.
Since an episode may end earlier due to collisions and successful goal
achievements, the number of RL steps per episode varies. Each training
consists of $200$ epochs with $30$ episodes per epoch, leading to a maximum of
$600\,000$ RL steps.

###  IV-B Randomized-goal experiment

The following discussion refers to the agent trained excluding and including a
safety shield as baseline and safe agent, respectively. We trained the
baseline agent with the same hyperparameters as the safe agent to achieve good
comparability. In order to improve the performance of the baseline agent, we
also tried to operate it at higher frequencies
($50\text{\,}\mathrm{H}\mathrm{z}$ and $250\text{\,}\mathrm{H}\mathrm{z}$) and
provide it with a high negative reward for collisions. Unfortunately, none of
the measures lead to successful training improvements, so we omit these
results. First, we want to highlight the results for completely randomized
goal positions. The collision rate comparison in Fig. 5a with and without a
safety shield clearly shows the advantage of our proposed method. The baseline
agent causes a safety-critical collision in approximately $20\text{\,}\%$ of
all episodes. Therefore, the baseline is too dangerous to be deployed in a
real-world environment, especially at the beginning of the training. In
contrast, the safety shield successfully prevents all safety-critical
collisions. Both agents show a similar success rate in Fig. 5b. Interestingly,
the baseline agent learns slightly faster than the safe agent despite its high
collision rate. At the end of the training, both agents have approximately the
same success rate of $65\text{\,}\%$. Fig. 5c evaluates how the remaining
$35\text{\,}\%$ of episodes end for the safe agent. Approximately
$10\text{\,}\%$ of episodes end in a safe collision caused by the static human
animation walking into the stationary robot. The remaining $25\text{\,}\%$ of
episodes ended by exceeding $T_{\text{max\\_episode}}$. A manual investigation
of the unsuccessful cases showed that most of the goals were not reachable
because humans blocked the path to the goal.

###  IV-C Human-evasion experiment

The human-evasion experiment highlights another advantage of our proposed
method. Since the human blocks the optimal path, the baseline agent collides
with the human in almost every episode in the first half of the training as
shown in Fig. 5d. In the second half of training, the baseline agent learned
to avoid collisions with the human more frequently, but still collides in
approximately $25\text{\,}\%$ of all episodes. Again, our safety shield
prevents all safety-critical collisions. Fig. 5e shows that the baseline agent
cannot reach the goal because it does not learn how to evade the human.
However, the safe agent learns to move around the human to reach the goal
effectively.

##  V CONCLUSIONS

Our results clearly show that our proposed safety shield is an effective
method to avoid collisions with humans in the working environment. Contrary to
existing RL methods for safe manipulator control, our safety shield provides
formal safety guarantees in highly dynamic and prior unknown human
environments. The human-evasion experiment also demonstrated that our safety
shield could significantly impact the training success in scenarios with a
high likelihood of a collision. We believe that our method enables the high-
level control of manipulators with RL agents in real-world human environments
while providing the necessary safety guarantees to any human operator. Our
subsequent goals are to refine the RL agent further and train it on a more
complex set of real-world human motions. Finally, we plan to test our safe RL
agent on actual robots and perform tasks in a human working environment.

## ACKNOWLEDGMENT

The authors gratefully acknowledge financial support by the Central Innovation
Programme of the German Federal Government under grant ZF4086004LP7 and the
Horizon 2020 EU Framework Project CONCERT under grant 101016007.

## References

  * [1] S. Gu, E. Holly, T. Lillicrap, and S. Levine, “Deep reinforcement learning for robotic manipulation with asynchronous off-policy updates,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2017, pp. 3389–3396.
  * [2] W. Gao, L. Graesser, K. Choromanski, X. Song, N. Lazic, P. Sanketi, V. Sindhwani, and N. Jaitly, “Robotic table tennis with model-free reinforcement learning,” in _Proc. of the IEEE/RSJ Int. Conf. on Intelligent Robots and Systems (IROS)_ , 2020, pp. 5556–5563.
  * [3] R. Li, A. Jabri, T. Darrell, and P. Agrawal, “Towards practical multi-object manipulation using relational reinforcement learning,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2020, pp. 4051–4058.
  * [4] E. Prianto, M. Kim, J.-H. Park, J.-H. Bae, and J.-S. Kim, “Path planning for multi-arm manipulators using deep reinforcement learning: Soft actor–critic with hindsight experience replay,” _Sensors_ , vol. 20, no. 20, pp. 1–22, 2020.
  * [5] M. El-Shamouty, X. Wu, S. Yang, M. Albus, and M. F. Huber, “Towards safe human-robot collaboration using deep reinforcement learning,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2020, pp. 4899–4905.
  * [6] T.-H. Pham, G. De Magistris, and R. Tachibana, “Optlayer - practical constrained optimization for deep reinforcement learning in the real world,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2018, pp. 6236–6243.
  * [7] B. Sangiovanni, G. P. Incremona, M. Piastra, and A. Ferrara, “Self-configuring robot path planning with obstacle avoidance via deep reinforcement learning,” _IEEE Control Systems Letters_ , vol. 5, no. 2, pp. 397–402, 2021\.
  * [8] “Robotics - safety requirements - part 1: Industrial robots,” International Organization for Standardization, Geneva, CH, Standard DIN EN ISO 10218-1:2021-09 DC, Sep. 2021.
  * [9] T. Haarnoja, A. Zhou, P. Abbeel, and S. Levine, “Soft actor-critic: Off-policy maximum entropy deep reinforcement learning with a stochastic actor,” in _Proc. of the 35th Int. Conf. on Machine Learning (ICML)_ , Jul. 2018, pp. 1861–1870.
  * [10] M. Andrychowicz, F. Wolski, A. Ray, J. Schneider, R. Fong, P. Welinder, B. McGrew, J. Tobin, P. Abbeel, and W. Zaremba, “Hindsight experience replay,” in _Proc. of the 31st Int. Conf. on Neural Information Processing Systems (NeurIPS)_ , 2017, pp. 5055–5065.
  * [11] Z. Ren, K. Dong, Y. Zhou, Q. Liu, and J. Peng, “Exploration via hindsight goal generation,” in _Proc. of the 33rd Int. Conf. on Neural Information Processing Systems (NeurIPS)_ , 2019.
  * [12] Z. Bing, M. Brucker, F. O. Morin, R. Li, X. Su, K. Huang, and A. Knoll, “Complex robotic manipulation via graph-based hindsight goal generation,” _IEEE Transactions on Neural Networks and Learning Systems_ , pp. 1–14, 2021, early access.
  * [13] J. Achiam, D. Held, A. Tamar, and P. Abbeel, “Constrained policy optimization,” in _Proc. of the 34th Int. Conf. on Machine Learning (ICML)_ , 2017, pp. 22–31.
  * [14] A. Ray, J. Achiam, and D. Amodei, “Benchmarking safe exploration in deep reinforcement learning,” _arXiv preprint arXiv:1910.01708_ , vol. 7, 2019\.
  * [15] G. Dalal, K. Dvijotham, M. Vecerik, T. Hester, C. Paduraru, and Y. Tassa, “Safe exploration in continuous action spaces,” _arXiv preprint arXiv:1801.08757_ , 2018.
  * [16] H. Krasowski, X. Wang, and M. Althoff, “Safe reinforcement learning for autonomous lane changing using set-based prediction,” in _Proc. of the IEEE 23rd Int. Conf. on Intelligent Transportation Systems (ITSC)_ , 2020, pp. 1–7.
  * [17] N. Hunt, N. Fulton, S. Magliacane, T. N. Hoang, S. Das, and A. Solar-Lezama, “Verifiably safe exploration for end-to-end reinforcement learning,” in _Proc. of the 24th Int. Conf. on Hybrid Systems: Computation and Control (HSCC)_ , 2021, article no. 14.
  * [18] Y. S. Shao, C. Chen, S. Kousik, and R. Vasudevan, “Reachability-based trajectory safeguard (RTS): A safe and fast reinforcement learning safety layer for continuous control,” _IEEE Robotics and Automation Letters_ , vol. 6, no. 2, pp. 3663–3670, 2021.
  * [19] N. Jansen, B. Könighofer, S. Junges, A. Serban, and R. Bloem, “Safe reinforcement learning using probabilistic shields,” in _Leibniz Int. Proc. in Informatics (LIPIcs)_ , vol. 171, 2020, pp. 31–316.
  * [20] M. Alshiekh, R. Bloem, R. Ehlers, B. Könighofer, S. Niekum, and U. Topcu, “Safe reinforcement learning via shielding,” in _Proc. of the 32nd AAAI Conf. on Artificial Intelligence (AAAI)_ , 2018.
  * [21] B. Könighofer, J. Rudolf, A. Palmisano, M. Tappler, and R. Bloem, “Online shielding for stochastic systems,” in _NASA Formal Methods_ , 2021, pp. 231–248.
  * [22] G. De Giacomo, L. Iocchi, M. Favorito, and F. Patrizi, “Foundations for restraining bolts: Reinforcement learning with LTLf/LDLf restraining specifications,” _Proc. of the 29th Int. Conf. on Automated Planning and Scheduling (ICAPS)_ , pp. 128–136, 2021.
  * [23] M. Althoff, A. Giusti, S. B. Liu, and A. Pereira, “Effortless creation of safe robots from modules through self-programming and self-verification,” _Science Robotics_ , vol. 4, no. 31, pp. 1–14, 2019.
  * [24] J. Achiam, “Spinning up in deep reinforcement learning,” OpenAI, 2018.
  * [25] B. D. Ziebart, “Modeling purposeful adaptive behavior with the principle of maximum causal entropy,” Ph.D. dissertation, Carnegie Mellon University, 2010\.
  * [26] D. Beckert, A. Pereira, and M. Althoff, “Online verification of multiple safety criteria for a robot trajectory,” in _Proc. of the IEEE 56th Conf. on Decision and Control (CDC)_ , 2017, pp. 6454–6461.
  * [27] A. Pereira and M. Althoff, “Calculating human reachable occupancy for guaranteed collision-free planning,” in _Proc. of the IEEE/RSJ Int. Conf. on Intelligent Robots and Systems (IROS)_ , 2017, pp. 4473–4480.
  * [28] ——, “Overapproximative human arm occupancy prediction for collision avoidance,” _IEEE Transactions on Automation Science and Engineering_ , vol. 15, no. 2, pp. 818–831, 2018.
  * [29] T. Kröger and F. M. Wahl, “Online trajectory generation: Basic concepts for instantaneous reactions to unforeseen events,” _IEEE Transactions on Robotics_ , vol. 26, no. 1, pp. 94–111, 2010.
  * [30] “Safety of machinery - positioning of safeguards with respect to the approach speeds of parts of the human body,” International Organization for Standardization, Geneva, CH, Standard DIN EN ISO 13855:2010-10 ST N, 2010.
  * [31] A. Pereira and M. Althoff, “Safety control of robots under computed torque control using reachable sets,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2015, pp. 331–338.
  * [32] S. Schepp, J. Thumm, S. B. Liu, and M. Althoff, “Sara: A tool for safe human–robot coexistence and collaboration through reachability analysis,” in _Proc. of the IEEE Int. Conf. on Robotics and Automation (ICRA)_ , 2022, early access.


# Asymptotic profiles of zero points of solutions to the heat equation

Hiroshi Ishii  Department of Mathematics, Hokkaido University, Sapporo,
060-0810, Japan<EMAIL_ADDRESS>
###### Abstract.

In this paper, we consider the asymptotic profiles of zero points for the
spatial variable of the solutions to the heat equation. By giving suitable
conditions for the initial data, we prove the existence of zero points by
extending the high-order asymptotic expansion theory for the heat equation.
This reveals a previously unknown asymptotic profile of zero points diverging
at $O(t)$. In a one-dimensional spatial case, we show the zero point’s second
and third-order asymptotic profiles in a general situation. We also analyze a
zero-level set in high-dimensional spaces and obtain results that extend the
results for the one-dimensional spatial case.

Date: . Corresponding author: H. Ishii

##  1\. Introduction

We consider a Cauchy problem

(1.1) $\begin{cases}\dfrac{\partial u}{\partial t}=\Delta_{x}u\quad(t>0,\
x\in{\mathbb{R}}^{d}),\vspace{3mm}\\\
u(0,x)=u_{0}(x)\quad(x\in{\mathbb{R}}^{d}),\end{cases}$

where $u=u(t,x)\in{\mathbb{R}}\ (t>0,x\in{\mathbb{R}}^{d})$, $\Delta_{x}$ is
Laplace operator for the spatial variable $x$, and $u_{0}\in
L^{1}({\mathbb{R}}^{d})\cap L^{\infty}({\mathbb{R}}^{d})$. Throughout this
paper, we assume that $\|u_{0}\|_{L^{1}}\neq 0$ and denote the Euclidean norm
and inner product on ${\mathbb{R}}^{d}$ by $|\cdot|$ and
$\left\langle\right.\cdot,\cdot\left.\right\rangle$, respectively.

In this paper, we consider the asymptotic behavior of the zero level set
${\mathcal{Z}}(t):=\\{x\in{\mathbb{R}}^{d}\mid u(t,x)=0\\}$, where $u(t,{x})$
is a unique bounded solution of (1.1) which can be written as follows:

(1.2)
$u(t,{x})=(G(t)*u_{0})(x)=\int_{{\mathbb{R}}^{d}}G(t,x-y)u_{0}(y)dy\quad(t>0,\
x\in{\mathbb{R}}^{d})$

of (1.1). Here, $G(t,{x})$ is the heat kernel defined as

(1.3) $G(t,{x}):=\dfrac{1}{(4\pi t)^{d/2}}e^{-|x|^{2}/4t}.$

The analysis of the zero level set is important to understanding the behavior
of sign-changing solutions of parabolic equations. To analyze the behavior of
sign-changing solutions, many researchers have focused on the zero level set
or the sign-changing number of the initial data [2, 4, 5, 14, 16]. Especially,
the case $d=1$ is particularly well analyzed. In parabolic equations, one of
the most fundamental properties is that ${\mathcal{Z}}(t)$ is a discrete set
[2]. Furthermore, it is known that the number of elements of
${\mathcal{Z}}(t)$ does not increase as time passes [5, 14].

In investigating the dynamics of the sign-changing solution in detail, the
asymptotic behavior of ${\mathcal{Z}}(t)$ has been analyzed to the heat
equation with $d=1$. Mizoguchi [16] estimated the upper bound of
${\mathcal{Z}}(t)$ in the case that $u_{0}(x)$ is sign-changing at finite
times and proved ${\mathcal{Z}}(t)\subset[-C_{0}t,C_{0}t]$ for sufficiently
large $t>0$ with some $C_{0}>0$. Moreover, it has been reported that there
exist $u_{0}(x)$ and $x^{*}(t)\in{\mathcal{Z}}(t)$ such that $x^{*}(t)>C_{1}t$
holds for sufficiently large $t>0$ with some $C_{1}>0$. Chung [4] analyzed the
asymptotic behavior of zero points when $u_{0}(x)$ belongs to
$L^{1}({\mathbb{R}},1+|x|^{k+1})$ for some $k\in{\mathbb{N}}$ and satisfies

(1.4) $\int_{{\mathbb{R}}}y^{j}u_{0}(y)dy=0\
(j=0,1,\ldots,k-1),\quad\int_{{\mathbb{R}}}y^{k}u_{0}(y)dy\neq 0.$

As the result of [4], for sufficiently large $t>0$, there exist
$x^{*}_{j}(t)\in{\mathcal{Z}}(t)\ (j=1,2,\ldots,k)$ such that

$\displaystyle\displaystyle\lim_{t\to+\infty}\dfrac{x^{*}_{j}(t)}{2\sqrt{t}}=h_{j},$

where $h_{j}\ (j=1,2,\ldots,k)$ are mutually different zero points of the
Hermite polynomial $H_{k}(x)$ with the order $k$ defined as

$\displaystyle H_{k}(x):=(-1)^{k}e^{x^{2}}\dfrac{d^{k}}{dx^{k}}[e^{-x^{2}}].$

The proof of these results is based on the Hermite polynomial approximation
deduced using the asymptotic expansion for the heat equation (see [4] and the
references therein). Chung further analyzed the behavior of the zero point of
the solution in the case $d\geq 2$ and also obtained results for the
spherically symmetric initial data.

We expect from these results that the elements of ${\mathcal{Z}}(t)$ have
various asymptotic profiles. However, even the case $d=1$ does not yield a
clear answer. Especially, even though there is a case that the asymptotic
behavior of an element of ${\mathcal{Z}}(t)$ is $O(t)$ as $t\to+\infty$, the
characterization of its coefficient has not been given. Indeed, it seems
interesting that the asymptotic behavior of some zero point is $O(t)$ while
$\sqrt{t}$ is the scale of self-similarity of the solution of the heat
equation. The aim of this paper is to explicitly give asymptotic profiles of
zero points and to reveal what properties of the initial data make a
difference in the orders.

As related problems, the level set of the solution to the heat equation has
been analyzed theoretically in connection with the hot spot problems [3, 9,
10, 12] and for application aspects such as level set methods [15] and image
processing [7, 8]. Among them, we mention the hot spot problem in the heat
equation. It is the study of the behavior of the maximum point of a solution,
and the maximum points are called hot spots. Under the non-negativity of the
initial condition and relatively general conditions, we can show that the hot
spots are critical points as pointed out in [3, 12]. Since the partial
derivatives of the solution are also solutions to the heat equation, analyzing
the behavior of the zero point also allows analysis of hot spots through
critical points. In this sense, the hot spot problem is an issue related to
the motivation of this paper. The hot spot problem for the heat equation on
unbounded domains in high dimensional space is often analyzed by imposing some
non-negativity on the initial data [3, 9, 10, 12]. To the best of our
knowledge, hot spots asymptotically move slower than $O(t)$ in many cases, and
there is no consideration of zero points moving at $O(t)$ which is the aim of
this paper.

In this paper, we develop higher-order asymptotic expansions of the heat
equation in order to analyze the behavior of the zero points in detail.
Furthermore, based on the analysis, the asymptotic behavior of the zero point
in the case $d=1$ is revealed up to the constant term. This result provides a
complete answer for the asymptotic behavior of the zero point in the one-
dimensional case when the initial data has fast decaying and finite sign
changes.

The paper is organized as follows: Section 2 presents results on higher-order
asymptotic expansions on moving coordinate systems for analyzing the
asymptotic behavior of zero points. In Section 3, the asymptotic behavior of
the zero point for $d=1$ is analyzed in detail with the proofs. We summarize
the analysis of specific cases for the higher dimensional case in Section 4.
Finally, we mention related problems in Section 5.

##  2\. Preliminaries

We prepare symbols and notations to explain the analytical methods used in
this study. For all multi index
$\alpha=(\alpha_{1},\ldots,\alpha_{d})\in({\mathbb{N}}\cup\\{0\\})^{d}$, we
define the sum of components as
$|\alpha|_{d}:=\displaystyle\sum^{d}_{j=1}\alpha_{j}$. Also, the $\alpha$
partial derivative to the smooth function $f(x)$ is defined as follows:

$\displaystyle\partial^{\alpha}_{x}f(x):=\dfrac{\partial^{|\alpha|}f}{\partial^{\alpha_{1}}_{x_{1}}\cdots\partial^{\alpha_{d}}_{x_{d}}}(x).$

We also define $x^{\alpha}=x^{\alpha_{1}}_{1}\cdots x^{\alpha_{d}}_{d}$ and
$\alpha!=\alpha_{1}!\cdots\alpha_{d}!$.

Under this preparation, it is well known that higher-order asymptotic
expansions for solutions such as the following:

######  Lemma 2.1.

[4, Theorem 2.2] Let $u$ be the bounded solution on
$(0,+\infty)\times{\mathbb{R}}^{d}$ to the heat equation (1.1) with initial
data $u_{0}\in L^{1}({\mathbb{R}}^{d};1+|x|^{n+1})$ for some nonnegative
integer $n$. There is then a positive constant $C=C(n,d)$ such that

$\displaystyle\left|u(t,x)-G(t,x)\displaystyle\sum_{|\alpha|\leq
n}\dfrac{\int_{{\mathbb{R}}^{d}}y^{\alpha}u_{0}(y)dy}{(4t)^{|\alpha|/2}\alpha!}\prod^{d}_{j=1}H_{\alpha_{j}}\left(\dfrac{x_{j}}{2\sqrt{t}}\right)\right|\leq
Ct^{-(n+d+1)/2}\||\cdot|^{n+1}u_{0}\|_{L^{1}}$

for all $x\in{\mathbb{R}}^{d}$.

We analyze the asymptotic behavior of the zero point by considering a natural
extension of this asymptotic expansion. The following conditions are given for
the initial data:

(2.1)
$\forall\lambda\in{\mathbb{R}}^{d},\quad\int_{{\mathbb{R}}^{d}}e^{\left\langle\right.\lambda,y\left.\right\rangle}|u_{0}(y)|dy<\infty.$

A concrete example satisfying (2.1) is bounded functions with compact support.

We assume that $u_{0}$ satisfies (2.1). For any $\eta\in{\mathbb{R}}^{d}$,
consider $v(t,x):=e^{<\eta,x>+t|\eta|^{2}}u(t,x+2t\eta)$. Then, $v(t,x)$ is a
solution of Cauchy problem

$\displaystyle\begin{cases}\dfrac{\partial v}{\partial t}=\Delta v\quad(t>0,\
x\in{\mathbb{R}}^{d}),\vspace{3mm}\\\
v(0,x)=e^{<\eta,x>}u_{0}(x)\quad(x\in{\mathbb{R}}^{d}).\end{cases}$

Since $v(0,x)$ belongs to $L^{1}({\mathbb{R}}^{d};1+|x|^{n+1})$ for all $n\geq
0$, we can apply Lemma 2.1. We define the bilateral Laplace transform as

$\displaystyle
U_{0}(\eta):=\int_{{\mathbb{R}}^{d}}e^{\left\langle\right.\eta,y\left.\right\rangle}u_{0}(y)dy,\quad\eta=(\eta_{1},\ldots,\eta_{d})\in{\mathbb{R}}^{d}.$

Notice that

$\displaystyle\partial^{\alpha}_{\eta}U_{0}(\eta)=\int_{{\mathbb{R}}^{d}}y^{\alpha}e^{\left\langle\right.\eta,y\left.\right\rangle}u_{0}(y)dy=\int_{{\mathbb{R}}^{d}}y^{\alpha}v(0,y)dy,$

we can deduce the following lemma.

######  Lemma 2.2.

Let $u$ be the bounded solution on $(0,+\infty)\times{\mathbb{R}}^{d}$ to the
heat equation (1.1) with initial data $u_{0}\in L^{1}({\mathbb{R}}^{d})\cap
L^{\infty}({\mathbb{R}}^{d})$ satisfying (2.1). For any $\eta\in{\mathbb{R}}$
and non-negative integer $n$, there is a positive constant $C=C(n,d)$ such
that $v(t,x)=e^{<\eta,x>+t|\eta|^{2}}u(t,x+2t\eta)$ satisfies

$\displaystyle\left|v(t,x)-G(t,x)\displaystyle\sum_{|\alpha|\leq
n}\dfrac{\partial^{\alpha}_{\eta}U_{0}(\eta)}{(4t)^{|\alpha|/2}\alpha!}\prod^{d}_{j=1}H_{\alpha_{j}}\left(\dfrac{x_{j}}{2\sqrt{t}}\right)\right|\leq
Ct^{-(n+d+1)/2}\|e^{\left\langle\right.\eta,\cdot\left.\right\rangle}|\cdot|^{n+1}u_{0}\|_{L^{1}}$

for all $x\in{\mathbb{R}}$.

Let us define the zero level set of $U_{0}$

$\displaystyle{\mathcal{N}}(U_{0}):=\\{\eta\in{\mathbb{R}}^{d}\mid
U_{0}(\eta)=0\\}.$

Assume that ${\mathcal{N}}(U_{0})\neq\emptyset$ and fix
$\eta^{*}\in{\mathcal{N}}(U_{0})$. Since $U_{0}(\eta)$ is a analytic function
on ${\mathbb{C}}^{d}$, there exists $k\in{\mathbb{N}}$ such that for all
$\alpha=(\alpha_{1},\ldots,\alpha_{d})\in({\mathbb{N}}\cup\\{0\\})^{d}$ with
$|\alpha|:=\displaystyle\sum^{d}_{j=1}\alpha_{j}<k$, we have
$\partial^{\alpha}_{\eta}U_{0}(\eta^{*})=0$, and there exists
$\tilde{\alpha}\in({\mathbb{N}}\cup\\{0\\})^{d}$ satisfying
$|\tilde{\alpha}|=k$ and $\partial^{\tilde{\alpha}}_{\eta}U_{0}(\eta^{*})\neq
0$. In this case, $\eta^{*}$ is called a zero point of $U_{0}(\eta)$ with
multiplicity $k$.

From Lemma 2.2, if ${\mathcal{N}}(U_{0})\neq\emptyset$ and
$\eta^{*}\in{\mathcal{N}}(U_{0})$ has multiplicity $k$, then there is a
$C=C(k,d,\eta^{*},u_{0})$ such that

$\displaystyle\left|v(t,x)-\dfrac{G(t,x)}{(4t)^{k/2}}\displaystyle\sum_{|\alpha|=k}\dfrac{\partial^{\alpha}_{\eta}U_{0}(\eta^{*})}{\alpha!}\prod^{d}_{j=1}H_{\alpha_{j}}\left(\dfrac{x_{j}}{2\sqrt{t}}\right)\right|\leq
Ct^{-(k+d+1)/2}$

for any $x\in{\mathbb{R}}^{d}$. This means that

(2.2)
$(4t)^{(d+k)/2}v(t,2\sqrt{t}x)\to\dfrac{e^{-|x|^{2}}}{\pi^{d/2}}\displaystyle\sum_{|\alpha|=k}\dfrac{\partial^{\alpha}_{\eta}U_{0}(\eta^{*})}{\alpha!}\prod^{d}_{j=1}H_{\alpha_{j}}\left(x\right)=:\tilde{H}_{d,k}(x;\eta^{*})$

holds uniformly as $t\to+\infty$. Since it is uniform convergence, we can
state the existence of zero points depending on the properties of
$\tilde{H}_{d,k}(x;\eta^{*})$. We have not yet found a unified method for this
issue, and thus the analysis must be done in a case-by-case manner.

######  Remark 2.3.

We note that the inequality in Lemma 2.2 can be rewritten as

$\displaystyle\left|v(t,x)-\displaystyle\sum^{n}_{m=1}\dfrac{1}{(4t)^{(d+m)/2}}\tilde{H}_{d,m}\left(\dfrac{x}{2\sqrt{t}};\eta^{*}\right)\right|\leq
Ct^{-(n+d+1)/2}\|e^{\left\langle\right.\eta,\cdot\left.\right\rangle}|\cdot|^{n+1}u_{0}\|_{L^{1}}.$

######  Remark 2.4.

The condition (2.1) can be discussed more weakly for the fast-decaying initial
data. However, it cannot be extended to the case of slow decay, as in
polynomial decay.

##  3\. The case $d=1$

Here we analyze the case $d=1$ in detail. We first consider the asymptotic
behavior of a zero point that is $O(t)$

######  Proposition 3.1.

Suppose that $u_{0}\in L^{1}({\mathbb{R}})\cap L^{\infty}({\mathbb{R}})$
satisfies (2.1). Assume that there exist $T>0$ and $x^{*}(t)\in Z(t)$ for
$t>T$ satisfying $x^{*}\in C(T,\infty)$ and

(3.1) $\limsup_{t\to+\infty}\left|\dfrac{x^{*}(t)}{2t}\right|<\infty.$

Then, $\theta:=\displaystyle\lim_{t\to+\infty}\dfrac{x^{*}(t)}{2t}$ exists and
belongs to ${\mathcal{N}}(U_{0})$.

###### Proof.

We first define

$\displaystyle\overline{\theta}:=\limsup_{t\to+\infty}\dfrac{x^{*}(t)}{2t},\quad\underline{\theta}:=\liminf_{t\to+\infty}\dfrac{x^{*}(t)}{2t}.$

and then fix $\theta_{0}\in[\underline{\theta},\overline{\theta}]$. Because
$x^{*}(t)\in C(T,\infty)$, there exists
$\\{t_{j}\\}_{j\in{\mathbb{N}}}\subset(T,\infty)$ such that $t_{j}\to+\infty$
and $x^{*}(t_{j})/2t_{j}\to\theta_{0}$ hold as $j\to+\infty$. From the
assumption that $x^{*}(t)\in{\mathcal{Z}}(t)$ for $t>T$, we have

$\displaystyle
u(t,x^{*}(t_{j}))=(G(t)*u_{0})(x^{*}(t_{j}))=0\quad\Leftrightarrow\quad\int_{{\mathbb{R}}}\exp\left(\frac{2x^{*}(t_{j})y-y^{2}}{4t_{j}}\right)u_{0}(y)dy=0$

for any $j\in{\mathbb{N}}$. From the dominated convergence theorem, we obtain

$\displaystyle
0=\int_{{\mathbb{R}}}e^{\theta_{0}y}u_{0}(y)dy=U_{0}(\theta_{0})$

as $j\to\infty$. This means that
$\theta_{0}=\overline{\theta}=\underline{\theta}$ and
$\theta_{0}\in{\mathcal{N}}(U_{0})$, because $U_{0}$ is an analytic function
on ${\mathbb{C}}$ and thus ${\mathcal{N}}(U_{0})$ is discrete. ∎

From this result, the asymptotic profiles corresponding to $O(t)$ are expected
to move at a certain speed without oscillating. In particular, if the upper
bound of ${\mathcal{Z}}(t)$ is $O(t)$, the asymptotic behavior of
${\mathcal{Z}}(t)$ is expected to change depending on the nature of
${\mathcal{N}}(U_{0})$.

We next mention the case that ${\mathcal{N}}(U_{0})=\emptyset$. Before
describing the result, we introduce a notation. For a function $f$ on
${\mathbb{R}}$ with $f(x)\not\equiv 0$, let $z(f)$ be the number of sign
changes; i.e. the supremum of $j$ such that

$\displaystyle f(x_{i})f(x_{i+1})<0,\quad i=1,2,\ldots,j$

for some $-\infty<x_{1}<x_{2}<\ldots<x_{j+1}<+\infty$.

If $f(x)$ satisfies $z(f)<+\infty$, then (3.1) in the statement of Proposition
3.1 always holds by Theorem 1.1 in [16]. We thus deduce the following result
from a proof by contradiction.

######  Corollary 3.2.

Suppose that $u_{0}\in L^{1}({\mathbb{R}})\cap L^{\infty}({\mathbb{R}})$
satisfies (2.1). Assume that ${\mathcal{N}}(U_{0})=\emptyset$ and
$z(u_{0})<\infty$. Then, there exists $T>0$ such that
${\mathcal{Z}}(t)=\emptyset$ for any $t>T$.

We analyze the case ${\mathcal{N}}(U_{0})\neq\emptyset$. By applying extended
asymptotic expansion, we obtain the following result:

######  Theorem 3.3.

Suppose that $u_{0}\in L^{1}({\mathbb{R}})\cap L^{\infty}({\mathbb{R}})$
satisfies (2.1). Assume that ${\mathcal{N}}(U_{0})\neq\emptyset$. Then, for
all $\eta^{*}\in{\mathcal{N}}(U_{0})$ with multiplicity $k\in{\mathbb{N}}$,
there exist $T>0$ and $x^{*}_{j}(t)\in{\mathcal{Z}}(t)\ (j=1,\ldots,k)$ for
$t>T$ satisfying

$\displaystyle
x^{*}_{j}(t)=2t\eta^{*}+2\sqrt{t}h_{j}+\dfrac{U^{(k+1)}_{0}(\eta^{*})}{(k+1)U^{(k)}_{0}(\eta^{*})}+o(1)\quad(t\to+\infty)$

for $j=1,\ldots,k$, where $h_{j}\ (j=1,\ldots,k)$ are mutually different zero
points of the Hermite polynomial $H_{k}(x)$ with the order $k$.

###### Proof.

Let $v(t,x)=e^{\eta^{*}x+t(\eta^{*})^{2}}u(t,x+2t\eta^{*})$, and analyze zero
points of $v(t,x)$ By applying Lemma 2.2 with $d=1$, $\eta=\eta^{*}$ and
$n=k+1$, there is $C>0$ such that

$\displaystyle\left|v(t,x)-G(t,x)\displaystyle\sum^{k+1}_{m=k}\dfrac{U^{(m)}_{0}(\eta^{*})}{(4t)^{m/2}m!}H_{m}\left(\dfrac{x}{2\sqrt{t}}\right)\right|\leq
Ct^{-(k+3)/2}$

for all $x\in{\mathbb{R}}$. Let $h$ be a zero point of $H_{k}(x)$. Then, we
rewrite the inequality as

$\displaystyle\left|v(t,2\sqrt{t}h+x)-G(t,2\sqrt{t}h+x)\displaystyle\sum^{k+1}_{m=k}\dfrac{U^{(m)}_{0}(\eta^{*})}{(4t)^{m/2}m!}H_{m}\left(h+\dfrac{x}{2\sqrt{t}}\right)\right|\leq
Ct^{-(k+3)/2}$

by simple computation and changing variable $x$ as $2\sqrt{t}h+x$. Note that
$C$ is re-taken here. Since we know that
$H^{\prime}_{k}(h)=2hH_{k}(h)-H_{k+1}(h)=-H_{k+1}(h)\neq 0$, we deduce

$\displaystyle(4t)^{(k+2)/2}G(t,2\sqrt{t}h+x)\displaystyle\sum^{k+1}_{m=k}\dfrac{U^{(m)}_{0}(\eta^{*})}{(4t)^{m/2}m!}H_{m}\left(h+\dfrac{x}{2\sqrt{t}}\right)$
$\displaystyle=\exp\left(-\dfrac{(x+2\sqrt{t}h)^{2}}{4t}\right)\left\\{2\sqrt{t}\dfrac{U^{(k)}_{0}(\eta^{*})}{k!}H_{k}\left(h+\dfrac{x}{2\sqrt{t}}\right)\right.$
$\displaystyle\hskip
170.71652pt\left.+\dfrac{U^{(k+1)}_{0}(\eta^{*})}{(k+1)!}H_{k+1}\left(h+\dfrac{x}{2\sqrt{t}}\right)\right\\}$
$\displaystyle\to
e^{-h^{2}}\left\\{\dfrac{U^{(k)}_{0}(\eta^{*})}{k!}H^{\prime}_{k}(h)x+\dfrac{U^{(k+1)}_{0}(\eta^{*})}{(k+1)!}H_{k+1}(h)\right\\}\quad(t\to+\infty)$
$\displaystyle=\dfrac{e^{-h^{2}}H_{k+1}(h)}{k!}\left\\{\dfrac{U^{(k+1)}_{0}(\eta^{*})}{k+1}-U^{(k)}_{0}(\eta^{*})x\right\\}.$

We note that the limit is right in the sense of local uniforms. We also know
the convergence rate of the limit, which allows us to construct a zero point
by using the intermediate value theorem. Thus, there exist $T>0$ and
$p^{*}(t)\in{\mathbb{R}}$ for $t>T$ such that

$\displaystyle v(t,p^{*}(t))=0,\quad
p^{*}(t)=2\sqrt{t}h+\dfrac{U^{(k+1)}_{0}(\eta^{*})}{(k+1)U^{(k)}_{0}(\eta^{*})}+o(1)\quad(t\to+\infty).$

Thus, there exist such zero point $p^{*}_{j}(t)$ of $v(t,x)$ for each
$h_{j}$,.

From the above discussion, setting $x^{*}_{j}(t)=2t\eta^{*}+p^{*}_{j}(t)\
(j=1,2,\ldots,k)$ yields $x^{*}_{j}(t)\in{\mathcal{Z}}(t)$ and the desired
asymptotic behavior. ∎

######  Remark 3.4.

We can compose $x^{*}_{j}(t)\in{\mathcal{Z}}(t)\ (j=1,\ldots,k)$ obtained in
Theorem 3.3 as $x^{*}_{j}\in C(T,\infty)\ (j=1,\ldots,k)$ by applying the
argument and results of [2].

From Theorem 3.3, for any $\eta^{*}\in{\mathcal{N}}(U_{0})$ and sufficiently
large $t>0$, there exists $x^{*}(t)\in{\mathcal{Z}}(t)$ such that

$\displaystyle\displaystyle\lim_{t\to+\infty}\dfrac{x^{*}(t)}{2t}=\eta^{*}.$

Thus, if $\eta^{*}\neq 0$, then $x^{*}(t)$ diverges with $O(t)$ as time
passes. In the case that $0\in{\mathcal{N}}(U_{0})$ with multiplicity $k>1$,
for sufficiently large $t>0$, there exists $x^{*}(t)\in{\mathcal{Z}}(t)$ such
that

$\displaystyle\displaystyle\lim_{t\to+\infty}\dfrac{x^{*}(t)}{2\sqrt{t}}=h,$

where $h$ is a zero point of $H_{k}(x)$. When $0\in{\mathcal{N}}(U_{0})$ with
multiplicity $k$ and $k$ is odd, for sufficiently large $t>0$, there is
$x^{*}(t)\in Z(t)$ satisfying

$\displaystyle\displaystyle\lim_{t\to+\infty}x^{*}(t)=\dfrac{U^{(k+1)}_{0}(0)}{(k+1)U^{(k)}_{0}(0)}=\dfrac{\displaystyle\int_{{\mathbb{R}}}y^{k+1}u_{0}(y)dy}{(k+1)\displaystyle\int_{{\mathbb{R}}}y^{k}u_{0}(y)dy}.$

Hence, by analyzing the zero points of $U_{0}(\eta)$ and their multiplicity,
we can understand the long-time behavior of elements of ${\mathcal{Z}}(t)$.

##  4\. The case $k=1$ in high dimensional space

The high-dimensional case is discussed. In this case, no unified results have
yet been obtained, and even the upper bound of the zero level set has not been
clarified. For the initial data $u_{0}(x)$ that is radial symmetric, Chung [4]
considered the case that $\eta^{*}=0\in{\mathcal{N}}(U_{0})$ and proved that
$\tilde{H}_{d,k}(x;\eta^{*})$ can be represented using the generalized
Laguerre polynomial. Chung also analyzed the case that
$\eta^{*}=0\in{\mathcal{N}}(U_{0})$ with mulitiplicity $k=1,2$ for $d=2$. This
reveals that there exists a zero point that is $O(\sqrt{t})$ even in higher
dimensions. However, when $\tilde{\eta}\neq 0$, the asymptotic behavior of the
zero point of the solution is not obtained. In particular, the asymptotic
profile of $O(t)$ is not known.

The analysis under general assumptions is a future problem, and here we
summarize the result of asymptotic profiles of zero points for the case $d\geq
2$, ${\mathcal{N}}(U_{0})\neq\emptyset$, and $k=1$, i.e. zero points of
$U_{0}(\eta)$ have multiplicity 1. Let $\eta^{*}\in{\mathcal{N}}(U_{0})$ with
multiplicity 1. Then, we obtain

$\displaystyle\tilde{H}_{d,1}(x,;\eta^{*})=\dfrac{2e^{-|x|^{2}}}{\pi^{d/2}}\displaystyle\sum^{d}_{j=1}\dfrac{\partial
U_{0}}{\partial\eta_{j}}(\eta^{*})x_{j}=\dfrac{2e^{-|x|^{2}}}{\pi^{d/2}}\left\langle\right.\nabla_{\eta}U_{0}(\eta^{*}),x\left.\right\rangle$

by using $H_{0}(p)=1,\ H_{1}(p)=2p$ for $p\in{\mathbb{R}}$. Then, the
following theorem is obtained as a result which is a natural extension of the
case $d=1$.

######  Theorem 4.1.

Suppose that $u_{0}\in L^{1}({\mathbb{R}}^{d})\cap
L^{\infty}({\mathbb{R}}^{d})$ satisfies (2.1). Assume that
${\mathcal{N}}(U_{0})\neq\emptyset$. Then, for all
$\eta^{*}\in{\mathcal{N}}(U_{0})$ with multiplicity $1$, there exist $T>0$ and
$x^{*}(t)\in{\mathcal{Z}}(t)$ for $t>T$ satisfying

$\displaystyle
x^{*}(t)=2t\eta^{*}+\dfrac{\Delta_{\eta}U_{0}(\eta^{*})}{2|\nabla_{\eta}U_{0}(\eta^{*})|^{2}}\nabla_{\eta}U_{0}(\eta^{*})+o(1)\quad(t\to+\infty).$

###### Proof.

By applying Lemma 2.2 with $n=2$ and multiplying $(4t)^{(d+2)/2}$ by both
sides of the inequality, we have

$\displaystyle\left|(4t)^{(d+2)/2}v(t,x)-\left\\{2\sqrt{t}\tilde{H}_{d,1}\left(\dfrac{x}{2\sqrt{t}};\eta^{*}\right)+\tilde{H}_{d,2}\left(\dfrac{x}{2\sqrt{t}}\right)\right\\}\right|\leq
Ct^{-1/2}$

for all $x\in{\mathbb{R}}^{d}$ with a positive constant $C>0$. Hence, for all
$x\in{\mathbb{R}}$,

(4.1)
$(4t)^{(d+2)/2}v(t,x)\to\left\langle\right.\nabla_{x}\tilde{H}_{d,1}\left(0;\eta^{*}\right),x\left.\right\rangle+H_{d,2}\left(0;\eta^{*}\right)$

holds as $t\to+\infty$ in the sense of uniform convergence. Since we deduce

$\displaystyle\dfrac{\partial\tilde{H}_{d,1}}{\partial
x_{m}}(x;\eta)=-\dfrac{4x_{m}e^{-|x|^{2}}}{\pi^{d/2}}\displaystyle\sum^{d}_{j=1}\dfrac{\partial
U_{0}(\eta^{*})}{\partial\eta_{j}}x_{j}+\dfrac{2e^{-|x|^{2}}}{\pi^{d/2}}\dfrac{\partial
U_{0}(\eta^{*})}{\partial\eta_{m}}$

for $m=1,2,\ldots,d$, we get
$\nabla_{x}\tilde{H}_{d,1}\left(0;\eta^{*}\right)=\dfrac{2}{\pi^{d/2}}\nabla_{\eta}U_{0}(\eta^{*})$.
Also, we have

$\displaystyle H_{d,2}(0)$ $\displaystyle=$
$\displaystyle\dfrac{1}{\pi^{d/2}}\displaystyle\sum_{|\alpha|=2}\dfrac{\partial^{\alpha}U_{0}(\eta^{*})}{\alpha!}\prod^{d}_{j=1}H_{\alpha_{j}}\left(0\right)$
$\displaystyle=$
$\displaystyle\dfrac{1}{2\pi^{d/2}}\displaystyle\sum^{d}_{j=1}\dfrac{\partial^{2}U_{0}(\eta^{*})}{\partial\eta^{2}_{j}}H_{2}\left(0\right)$
$\displaystyle=$
$\displaystyle-\dfrac{1}{\pi^{d/2}}\Delta_{\eta}U_{0}(\eta^{*}).$

Hence, the convergence of (4.1) is rewritten as

$\displaystyle\left\langle\right.\nabla_{x}\tilde{H}_{d,1}\left(0;\eta^{*}\right),x\left.\right\rangle+H_{d,2}\left(0;\eta^{*}\right)=\dfrac{1}{\pi^{d/2}}\\{2\left\langle\right.\nabla_{\eta}U_{0}(\eta^{*}),x\left.\right\rangle-\Delta_{\eta}U_{0}(\eta^{*})\\}.$

Then,
$x=\dfrac{\Delta_{\eta}U_{0}(\eta^{*})}{2|\nabla_{\eta}U_{0}(\eta^{*})|^{2}}\nabla_{\eta}U_{0}(\eta^{*})\in{\mathbb{R}}^{d}$
satisfies

$\displaystyle
2\left\langle\right.\nabla_{\eta}U_{0}(\eta^{*}),x\left.\right\rangle-\Delta_{\eta}U_{0}(\eta^{*})=0.$

Since (4.1) is uniformly convergent and the convergence rate is known, we can
construct the zero point by using the intermediate value theorem and the
implicit function theorem. Therefore, there exist $T>0$ and
$p^{*}(t)\in{\mathcal{Z}}(t)$ for $t>T$ satisfying

$\displaystyle v(t,p^{*}(t))=0\ (t>T),\quad
p^{*}(t)=\dfrac{\Delta_{\eta}U_{0}(\eta^{*})}{2|\nabla_{\eta}U_{0}(\eta^{*})|^{2}}\nabla_{\eta}U_{0}(\eta^{*})+o(1)\quad(t\to+\infty).$

Thus, setting $x^{*}(t)=2t\eta^{*}+p^{*}(t)$ yields
$x^{*}(t)\in{\mathcal{Z}}(t)$ and the desired asymptotic behavior. ∎

Here we treat the case of the radial symmetric initial data as a corollary.
Let $u_{0}$ be radial symmetric function on ${\mathbb{R}}^{d}$ satisfying
(2.1). Then, $U_{0}(\eta)$ is also a radial symmetric function on
${\mathbb{R}}^{d}$. We assume that ${\mathcal{N}}(U_{0})\neq\emptyset$. Since
$U_{0}(\eta)$ is an analytic function, ${\mathcal{N}}(U_{0})$ must be discrete
with respect to the radial direction. Then, $U_{0}(\eta)=U_{s}(|\eta|)$
satisfies

$\displaystyle\dfrac{\partial
U_{0}}{\partial\eta_{j}}(\eta)=\dfrac{\eta_{j}}{|\eta|}U^{\prime}_{s}(|\eta|),\quad\dfrac{\partial^{2}U_{0}}{\partial\eta^{2}_{j}}(\eta)=\dfrac{\eta^{2}_{j}}{|\eta|^{2}}U^{\prime\prime}_{s}(|\eta|)+\left(\dfrac{1}{|\eta|}-\dfrac{\eta^{2}_{j}}{|\eta|^{3/2}}\right)U^{\prime}_{s}(|\eta|)$

for $|\eta|>0$ and $j=1,2,\ldots,d$. Hence,

$\displaystyle\dfrac{\Delta_{\eta}U_{0}(\eta)}{2|\nabla_{\eta}U_{0}(\eta)|^{2}}\nabla_{\eta}U_{0}(\eta)=\dfrac{1}{2}\left(\dfrac{U^{\prime\prime}_{s}(|\eta|)}{U^{\prime}_{s}(|\eta|)}+\dfrac{d-1}{|\eta|}\right)\dfrac{\eta}{|\eta|}$

holds if $U^{\prime}_{s}(|\eta|)\neq 0$. This yields the following result:

######  Corollary 4.2.

Suppose that $u_{0}\in L^{1}({\mathbb{R}}^{d})\cap
L^{\infty}({\mathbb{R}}^{d})$ is a radial function satisfying (2.1). Assume
that ${\mathcal{N}}(U_{0})\neq\emptyset$. Then, for any
$\eta^{*}\in{\mathcal{N}}(U_{0})$ with multiplicity $1$, there exist $T>0$ and
$r^{*}(t)>0$ for $t>T$ satisfying

$\displaystyle\\{x\in{\mathbb{R}}^{d}\mid|x|=r^{*}(t)\\}\subset{\mathcal{Z}}(t)$

for $t>T$ and

$\displaystyle
r^{*}(t)=2t|\eta^{*}|+\dfrac{1}{2}\left(\dfrac{U^{\prime\prime}_{s}(|\eta^{*}|)}{U^{\prime}_{s}(|\eta^{*}|)}+\dfrac{d-1}{|\eta^{*}|}\right)+o(1)\quad(t\to+\infty).$

######  Remark 4.3.

Note that $k>1$ when $\eta^{*}=0\in{\mathcal{N}}(U_{0})$. In such a case,
Chung [4] has obtained a result on the behavior of zero points moving in
$O(\sqrt{t})$.

When $k\geq 2$, we may see the behavior of term corresponding to $O(\sqrt{t})$
depending on the nature of $\tilde{H}_{d,k}(x;\eta^{*})$, but it is difficult
to find general properties about the zero point. There may be some good ways
to analyze them, but we have not found them yet. This point is currently left
as an open problem.

##  5\. Related problems

Finally, some related problems are mentioned as developmental questions.

###  5.1. Hot spot problems

In Euclidean space ${\mathbb{R}}^{d}$, it is known that when the initial data
is non-zero and non-negative, the hot spot converges to a point determined by
the initial data. When the initial data is sign-changing, the asymptotic
behavior of the hot spot is different from the case that the initial data is
non-negative. In fact, we can construct an example where all diverge at the
zero point of $u_{x}(t,x)$, which is a solution of the heat equation in the
case of $d=1$. This means that it corresponds to the case where all the
critical points, which include the maximum and minimum points of the solution
of the heat equation, diverge. To reveal the motion of hot spots, a more
detailed evaluation is needed to determine which critical points will become
hot spots. The hot spot problem for the initial data that is sign-changing is
one interesting future problems.

###  5.2. Reaction-diffusion equations

It is also interesting to investigate the dynamics near the constant state of
the reaction-diffusion equation:

$\displaystyle\dfrac{\partial u}{\partial t}=\Delta_{x}u+f(u)\quad(t>0,\
x\in{\mathbb{R}}^{d}),$

where $u=u(t,x)\in{\mathbb{R}}$ and $f:{\mathbb{R}}\to{\mathbb{R}}$ is a
smooth function. If $f(u)$ satisfies $f(0)=0$, then we can expect that the
solution dynamics can be approximated by a solution of a linearized equation
around $u=0$ as long as the solution is sufficiently small. Hence, the
behavior of the zero point of the solution can be analyzed by attributing it
to the diffusion equation, and the influence of the nonlinear term may be
seen.

One particularly interesting topics is the case of Fujita-type nonlinear terms
$f(u)=|u|^{p-1}u$ for $p>1$. Mizoguchi [16] discussed the upper bound of the
zero level set when $d=1$. It reported that the zero level set has similar
properties to the diffusion equation when $p>3$ and the decay of the global
solution is sufficiently fast. Furthermore, Mizoguchi and Yanagida [17]
reported that a global solution exists for $d=1$ by changing the sign of the
initial data even when $p$ is smaller than the Fujita exponent. They proposed
that a sign change weakens the condition for the global existence of the
solution to the Fujita equation. The analysis of the zero level set may
therefore be an indicator for considering the conditions for the global
existence of a solution. Analysis using higher-order asymptotic expansions of
solutions such as this paper may be useful and should be considered in the
future.

###  5.3. Nonlocal diffusion equations

The asymptotic behavior of the zero-level set of another diffusion equation
can be considered as another extension of the problem. As an example, we
consider the case of the space fractional diffusion equation:

$\displaystyle\begin{cases}\dfrac{\partial u}{\partial
t}+(-\Delta_{x})^{s}u=0\quad(t>0,\ x\in{\mathbb{R}}^{d}),\vspace{3mm}\\\
u(0,x)=u_{0}(x)\quad(x\in{\mathbb{R}}^{d}),\end{cases}$

where $(-\Delta)^{s}$ is the fractional power of the Laplace operator with
$0<s<1$ In this case, it can be deduced that the zeros disappear at a finite
time when $u_{0}\in L^{1}({\mathbb{R}}^{d})$ is compactly supported and its
integral value is non-zero, because the evaluation of the relative error has
been derived [13]. In addition, when the integral value of the initial data is
zero, the asymptotic behavior of the zero level set can be considered. The
specific behavior of the zero level set can be considered using the asymptotic
expansion of the solution [11].

Furthermore, the asymptotic behavior of the zero point of the solutions to the
time-fractional diffusion equation (the diffusion-wave equation) [6] and the
nonlocal diffusion equation [1] may also be considered, but we leave them as
open problems.

## Acknowledgments

The author expresses his sincere gratitude to Prof. Shin-Ichiro Ei (Hokkaido
University), Prof. Eiji Yanagida (Tokyo Institute of Technology) and Prof.
Kazuhiro Ishige (The University of Tokyo) for their useful suggestions and
valuable advice and thanks to Edanz (https://jp.edanz.com/ac) for editing a
draft of this manuscript. The author is partially supported by JSPS KAKENHI
Grant Number JP21J10036 and JP23K13013.

## References

  * [1] F. Andreu-Vaillo, J. Mazón, J. Rossi, J. Toledo-Melero,  Nonlocal diffusion problems,  Math. Surveys Monogr. 165, AMS, Providence, RI, 2010.
  * [2] S. Angenent,  The zero set of a solution of a parabolic equation,  J. Reine Angew. Math. 390 (1988), 79–96.
  * [3] I. Chavel and L.Karp,  Movement of hot spots in Riemannian manifolds,  J. Analys. Math., 55 (1990), 271-286.
  * [4] J. Chung,  Long-time asymptotics of the zero level set for the heat equation,  Quarterly of Applied Mathematics, 70(4) (2012) 705–720.
  * [5] A. Ducrot, T. Giletti, H. Matano Existence and convergence to a propagating terrace in one-dimensional reaction-diffusion equations Trans. Amer. Math. Soc., 366 (2014), 5541-5566.
  * [6] S. D. Eidelman, A. N. Kochubei,  Cauchy problem for fractional diffusion equations,  J. Differential Equations, 199(2) (2004) 211–255.
  * [7] M.-H. Giga, Y. Giga, T. Ohtsuka, N. Umeda,  On behavior of signs for the heat equation and a diffusion method for data separation,  Communications on Pure and Applied Analysis, 12(5) (2013) 2277-2296.
  * [8] T. Hidaka, K. Imamura, T. Hioki et al.,  Prediction of compound bioactivities using heat-diffusion equation,  Patterns, 1(9) (2020) 100140.
  * [9] K. Ishige,  Movement of hot spots on the exterior domain of a ball underthe Neumann boundary condition,  J. Differential Equations, 212(2) (2005) 394–431.
  * [10] K. Ishige,  Movement of hot spots on the exterior domain of a ball under the Dirichlet boundary condition,  Adv. Differential Equations, 12(10) (2007) 1135–1166.
  * [11] K. Ishige, T. Kawakami, H. Michihisa,  Asymptotic expansions of solutions of fractional diffusion equations,  SIAM J. Math. Anal. 49(3) (2017) 2167–2190.
  * [12] S. Jimbo, S. Sakaguchi,  Movement of hot spots over unbounded domains in ${\mathbb{R}}^{N}$,  Journal of Mathematical analysis and Applications, 182 (1994) 810–835.
  * [13] J. L. Vázquez,  Asymptotic behaviour for the fractional heat equation in the Euclidean space,  Complex Var. Elliptic Equ. 63(7-8) (2018) 1216–1231.
  * [14] H. Matano,  Nonincrease of the lap-number of a solution for a one-dimensional semilinear parabolic equation,  J. Fac. Sci. Univ. Tokyo Sect.1A Math. 29(2) (1982) 401–441.
  * [15] B. Merriman, J. Bence, S. Osher Motion of Multiple Junctions: A Level Set Approach Journal of Computational Physics 112 (1994) 334–363.
  * [16] N. Mizoguchi,  Asymptotic behavior of zeros of solutions for parabolic equations,  Journal of Differential Equations, 170(1) (2001) 51–67.
  * [17] N. Mizoguchi, E. Yanagida,  Critical exponents for the blow-up of solutions with sign changes in a semilinear parabolic equation,  Math. Ann., 307 (1997) 663–675.