Citation: Pandadaran, S.H., H. Yanagisawa, B. Shibazaki, Y. Fujii (2022), Solving the puzzle of the 1996
Biak, Indonesia tsunami, Synopsis of IISEE-GRIPS Master’s Thesis.

SOLVING THE PUZZLE OF THE 1996 BIAK,
INDONESIA TSUNAMI
Sidiq Hargo Pandadaran1
MEE21719

Supervisor: Hideaki YANAGISAWA2,
Bunichiro SHIBAZAKI3, Yushiro FUJII3
ABSTRACTOn February 17, 1996, an earthquake occurred northeast of Biak Island, Indonesia, and caused a tsunami.
Interestingly, the southwest side of Biak Island, which was not facing the epicenter, had a higher tsunami
runup than the facing side (Imamura et al., 1997). Previous researchers assumed that the earthquake
triggered submarine landslides. As no one has addressed this phenomenon, this assumption remains an
unsolved hypothesis.
The tsunami arrival time obtained from local people’s eyewitness testimonies was used toperform backward tsunami raytracing. Considering the raytracing result and multibeam topography, we
found two submarine landslide candidates: one large submarine landslide and a small submarine
landslide were located in the southeast and south of Biak Island, respectively. Since the small submarine
landslide only had a small effect on the land, we only performed a tsunami simulation for the large
submarine landslide. The result showed that the submarine landslide located at 135.62oE and -1.01oS
with a geometry of about 950 m x 5000 m and a thickness of about 50 m seems to explain the observed
runup and arrival time.
Previous researchers made a slip distribution without considering the submarine landslide
event. As a result, their model could not explain the observed runup in the southwest coastal area of Biak
Island. To accommodate this problem, we propose a new model by combining the submarine landslide
model with a modified fault model from the previous researchers. Our new model explains observed
runup heights well; we obtained a geometric mean of 1.00 and a geometric standard deviation of 1.40.
Keywords: Tsunami simulation, submarine landslide, backward tsunami raytracing, Biak.
1. INTRODUCTION
On 17 February 1996 at 05:59 UTC, an unpredicted earthquake with Mw 8.2 occurred in the northeast
of Biak Island and caused a tsunami (Figure 1). Interestingly, an earthquake with a thrust fault
mechanism occurred in the NGT segment, where there was no historical record of a large earthquake inthis area (Okal, 1999). The tsunami generated not only hit Biak Island but also spread to Taiwan, Japan,
Canada, and the USA.
One month after the event, the International Tsunami Survey Team (ITST) did a tsunami
survey on Biak Island (Imamura et al., 1996; Imamura et al., 1997). The tsunami hit the entire coast of
Biak Island, with an average height of 2 m. However, surprisingly, a tsunami as high as 7.7 m occurredin Mardori, the southwest part of Biak Island. This location did not face the earthquake’s epicenter but
had a tsunami runup higher than that on the northeast coast, with a maximum runup of 5.4 m in Korem(Figure 1).

1

Indonesian Agency for Meteorology, Climatology, and Geophysics (BMKG), Indonesia.
Associate Professor, Tohoku Gakuin University.
3
International Institute of Seismology and Earthquake Engineering, Building Research Institute.
2

1

Earthquake

Submarine Landslide

Figure 1. Map of the study area.

Model
Imamura et al. (1997)
Tanioka and Okada (1997)
Henry and Das (2002)
Arimuko et al. (2020)

Table 1. Fault model used in this study.
Length Width Strike Dip Rake
(km)
(km)
(o)
(o)
(o)
180
50
112
11
91
200
75
115
10
90
320
100
109
9
72
272
110
103
11
69

Slip
(m)
8
3
4
0.922

Mag.
7.9
8.1
8.2
8.05

Depth
(km)
15
10
9.5
13.7

Several researchers have conducted research on tsunami surveys, tsunami damage, and runup
simulation for Biak Island (Imamura et al., 1996; Imamura et al., 1997; Matsutomi et al., 2001), but
none of them have solved the unexpected highest runup at southwest coastal area. Until now, thesubmarine landslide that contributed to the 1996 Biak, Indonesia tsunami is still only an assumption.
Therefore, through this research, we want to solve this hypothesis.
2. DATA
2.1. Bathymetry and Topography
A three-level nested grid layer that consists of 10-arcsecond resolutions for each layer was used in thisstudy (Figure 1). Batimetri Nasional (BATNAS) is used in layer 1. A nautical chart, high resolution
multibeam, and Digital Elevation Model Nasional (DEMNAS) are combined for the second and thirdlayer nested grids.
2.2. Tsunami Travel Time, Runup, and Tectonic Source
Due to the absence of tide gauge stations around Biak Island, Matsutomi et al. (2001) made the estimatedtsunami arrival time by interviewing local people. Runup height (Imamura et al., 1996) was collected

2

(a)

Northeast(b)South-Southwest(c)South-Southwest(d)Northeast(e)South-Southwest

Figure 2. Comparison of observed and calculated tsunami runups. (a) and (b) are from fault models.
(c) is from candidates 1 and 2 of the submarine landslide sources. (d) and (e) are from combination
between candidate 1 of the submarine landslide and fault model modified from Imamura et al. (1997).
during the field survey (Table 1 in Imamura et al. (1996)). On the northeast coastal area, facing theepicenter, the highest runup with the 5.4 m height was recorded at Korem. On the other hand, the tsunami
runup significantly increased in the backside area from the epicenter, where the highest tsunami runup
recorded was 7.7 m at Mardori.
To identify the possibility of submarine landslide, we applied the same technique as Okal etal. (2002). We fitted the polynomial function for the northeast and south-southwest coastal areas byusing the “trend1d” command in Generic Mapping Tools 6 (GMT; Wessel et al., 2019) and smoothing
manually (Figure 1), then compared it with Figure 2 in Okal et al. (2002). From this result, we canconclude that the northeast and southwest coastal areas were affected by the tsunami from earthquakeand submarine landslides, respectively.
The fault models from previous researchers (Table 1) were used in the tsunami simulation todetermine whether the fault model can explain the distribution of the tsunami runup or not. Since there
was no information for the detailed subfault parameter from Henry and Das (2002) and Arimuko et al.
(2020), for simplicity, we assumed their model as a single fault. Also, we set the depth of 10 km for the
Tanioka and Okada’s fault (1997) because they did not mention it in their paper.

3

(a)

(b)

Figure 3. (a) Possible submarine landslide area in the south-southwest of Biak Island from backward
tsunami raytracing and (b) detailed geomorphology near the southwest possible area.
3. METHODOLOGY
3.1. Backward Tsunami Raytracing
We performed backward tsunami raytracing to find out the location of the possible submarine landslidein the southwest coastal area of Biak Island. The sea-level station in the coastal area was used to draw
the inverse arrival timeline based on the tsunami arrival from the local people's eyewitnesses (Table 2
in Matsutomi et al. (2001)). Since there was no accurate location where the people saw the tsunami
coming, we assumed that the observation point was near the coordinates of the tsunami runup survey(Figure 1). For the travel time estimation, we used Tsunami Travel Time (TTT) software package (TTTSDK v 4.0.1) provided by ITIC (2021).
3.2. Tsunami Simulation
The tsunami wave generation in the southwest coastal area is assumed to be a submarine landslide. To
estimate the initial sea surface affected by sliding mass, we adopted a method that is similar to a slump(Watts et al., 2005). For tsunami propagation, we used TUNAMI compiled by Yanagisawa (2022). The
nonlinear longwave equation on cartesian coordinate is solved using staggered leapfrog finite-difference
scheme with the dispersive term.
3.3. Result Validation
K and κ index by Aida (1978) have been commonly used to validate the tsunami simulation. This
technique uses geometric mean, K, and geometric standard deviation, κ, to assess the result. The Japan
Society of Civil Engineers (JSCE, 2002) made a threshold for a valid tsunami source model with thevalue of 0.95 < K < 1.05 and κ < 1.45.
4. RESULTS AND DISCUSSION
4.1. Tsunami Simulation from Fault Model
Two-layer nested grids with 10 arcseconds were applied in the tsunami simulation. We could not use amore detailed spatial grid size since we found that some of the observation data coordinates wereincorrect because they were located on the shallow sea area, not the land. For the result, all fault models
could generally explain the distribution of the tsunami runup on the northeast coast (Figure 2a). However,none of the models could explain the 7.7-meter tsunami runup on the southwest coast (Figure 2b). Thus,
these results strengthen the hypothesis regarding the existence of submarine landslides.

4

(a)

(b)

Figure 4. Combination model between submarine landslide and fault model modified from Imamura etal. (1997). (a) Initial condition and (b) maximum tsunami height.
4.2. Backward Tsunami Raytracing and Submarine Landslide Location
We calculated the backward tsunami raytracing using TTT and 6 grid arcseconds of constructedbathymetry (Figure 3a). To find out how many submarine landslides contributed to this event, we
compared the distribution of tsunami runup in the south-southwest coastal area (Figure 1) against the2018 Palu, Indonesia tsunami case (Figure 1 in Nakata et al. (2020)). Most submarine landslides
occurred near the east coast of Palu Bay and produced many peaks in the runup distribution. However,in our case, only one peak was clearly visible near Mardori (Figure 1). This means that only one
submarine landslide contributed to this event.
In Addition, the shorter arrival time in the south part of Biak Island indicated another possible
submarine landslide (Figure 3a). Taking the width of the sea floor collapse, the runup distribution, andthe highest runup was located in Sorido; therefore, we concluded that the submarine landslide is small
and facing the Sorido (Sorido will be the backside area of this submarine landslide).
We analyzed and delineated the seafloor's topography to identify the steep cliff and the contour
profile related to the submarine landslide (Figure 3b). From the shallow to the deep sea floor, there areterraced fault scarps, and there is a canyon gully. Because the available multibeam data only covered
half of the possible area, we suspected that there were two candidates for the large submarine landslide,which were located near the steep cliff in the left and right of the canyon gully with the locations of
135.62oE-1.01oS and 135.66oE-1.08oS (Figure 3b), respectively.
4.3. Tsunami Simulation from Submarine Landslide and Fault Model
The submarine landslide parameters were determined based on the slope geometry and trial and error.
The first candidate was easier to be identified, given the slight disturbance in the contour below theslope, which indicates a submarine landslide deposit. Since there was no evident proof in the second
candidate, for simplicity, we assumed that the parameter for the second candidate was the same as thatof the first candidate.
The submarine landslide tsunami simulation was only conducted in the southwest side source
because for the south of Biak Island case, the submarine landslide was small, did not significantly impact
the land, and the runup distribution could be described well by the fault model (Figure 2b). Thus, it
could be a bias if we simulated the submarine landslide model because we did not know whether therecorded runups were produced by the tsunami from the earthquake or a small submarine landslide.
The first candidate got better K and κ values of 1.00 and 1.54 than the second candidate, with
only 0.95 and 1.64 with both submarine landslide geometry of 950 m x 5000 m and a thickness of 50 m.
From the result above, we conclude that the most possible submarine landslide source is the firstcandidate. To create a model that could describe the runup distribution on all sides of Biak Island, we
did a tsunami simulation by combining the fault model (Table 1) and the submarine landslide parameter5

of the first candidate. The results show that the model of Imamura et al. (1997) has the best results
among the others with the K and κ values of 0.92 and 1.41, respectively.
However, with a K value of 0.92, it shows that this combination of models was stilloverestimated and did not meet the requirements of the valid tsunami source from JSCE (2002). To
accommodate that problem, we simplify to reduce the amount of slip in their model. Finally, we got theideal slip of 7.1 m with the K and κ of 1.00 and 1.40, and the runup distribution from our new modeland the detailed tsunami simulation model can be seen in Figures 2a, 2e, 4a, and 4b.
5. CONCLUSION
We agree with the hypothesis of previous researchers that the 1996 Biak, Indonesia earthquake triggered
several submarine landslides surrounding Biak Island. Based on our analysis, this earthquake triggered
one big submarine landslide and a small submarine landslide in the southwest and south of Biak Island,
respectively. We used a submarine landslide model (Watts et al., 2005) to calculate the tsunami runup
caused by the bigger one because the small submarine landslide did not impact the land significantly.
Our simulation found that the submarine landslide located in the 135.62oE and -1.01oS with thedimension of about 950 m x 5000 m and a thickness of about 50 m seems to be able to explain theobserved runup and arrival time. We proposed a new model by combining and modifying the submarine
landslide parameter and the previous researcher's model. We chose Imamura et al. (1997) model because
their model had good K and κ values among the other models. To reproduce the best fit model, we
modified the amount of slip in their model from 8 m to 7.1 m, and we got K 1.00 and κ 1.40 for our newmodel.
ACKNOWLEDGEMENTS
I want to acknowledge my sincere gratitude to my supervisors, Dr. Hideaki Yanagisawa, Dr. Bunichiro
Shibazaki, and Dr. Yushiro Fujii, for all the guidance and advice through this work and also to all otherlecturers and staffs at Building Research Institute (BRI), Japan International Corporation Agency (JICA),
and National Graduate Institute for Public Policy (GRIPS) for their help and support during my study
in Japan.
REFERENCES
Aida, I., 1978, Journal of Physics of the Earth, 26, 57-73.
Arimuko, A., et al., 2020, Journal of Physics: Conference Series, 1805, 012017.
Henry, C., and Das, S., 2002, Journal of Geophysical Research: Solid Earth, 107, 11.
Imamura, F., et al., 1996, Eos: Transactions American Geophysical Union, 78, 197-201.
Imamura, F., et al., 1997, Report on the survey of Irian Jaya earthquake and tsunami, Indonesia of 17
February 1996 by International Tsunami Survey Team.
Japan Society of Civil Engineers, 2002, Tsunami assessment method for nuclear power plants inJapan, http://www. jsce.or.jp/committee/ceofnp/Tsunami/eng/JSCE_Tsunami_060519.pdf
Matsutomi, H., et al., Natural Hazards, 24, 199-212.
Nakata, K., et al., 2020, Earth, Planets and Space, 72, 1-16.
Okal, E.A., 1999, Seismogenic and Tsunamigenic Processes in Shallow Subduction Zones, 633-675.
Tanioka, Y., and Okada, M., 1997, The International Journal of The Tsunami Society, 15, 67-79.
Watts, P., et al., 2005, Journal of Waterway, Port, Coastal, and Ocean Engineering, 131, 298-310.
International Tsunami Information Center, 2021, Tsunami Travel Time Software, http://itic.iocunesco.org/index.php?option=com_content&view=article&id=1576&Itemid=2933.
Wessel, P., et al., 2019, Geochemistry, Geophysics, and Geosystems, 20, 5556–5564.
Yanagisawa, H., 2021-2022, Lecture notes on Numerical Simulation Tsunami Modelling, IISEE/BRI.

6