# Laser slicing: a thin film lift-off method for GaN-on-GaN technology

Vladislav Voronenkov1,2<EMAIL_ADDRESS>Natalia Bochkareva2 Ruslan
Gorbunov1,2 Andrey Zubrilov1,2 Viktor Kogotkov1 Philipp Latyshev1 Yuri
Lelikov1,2 Andrey Leonidov1 Yuri Shreter1,2<EMAIL_ADDRESS>1TRINITRI-
Technology LLC, St Petersburg 194223, Russia
2Ioffe Institute, St Petersburg 194021, Russia

###### Abstract

A femtosecond laser focused inside bulk GaN was used to slice a thin GaN film
with an epitaxial device structure from a bulk GaN substrate. The demonstrated
laser slicing lift-off process did not require any special release layers in
the epitaxial structure. GaN film with a thickness of 5 $\mu$m and an InGaN
LED epitaxial device structure was lifted off a GaN substrate and transferred
onto a copper substrate. The electroluminescence of the LED chip after the
laser slicing lift-off was demonstrated.

laser slicing; lift-off; InGaN LED; GaN-on-GaN

Bulk GaN substrate is a crucial element in the production of reliable, high-
current density GaN-based devices. However, the widespread use of GaN
substrates is currently limited due to their high cost. Separation of the
epitaxial device structure from the bulk GaN substrate allows the reuse of the
expensive substrate multiple times and thus reduces the effective cost. Also,
the separation of the epitaxial structure is necessary to create
heterogeneously integrated devices heterogeneous-gan-cmos-2009 ;
HETEROGENEOUS-LED-2015 ; HETEROGENEOUS-2010 , as well as to improve heat
dissipation DIAMOND-HEMT-2014-Kuball ; HETEROGENEOUS-LED-2015  and light
extraction LLO-InGaN-AKASAKI2014 .

To lift-off GaN films from foreign substrates multiple methods were developed:
laser lift-off LLO-Kelly-1996 ; natural stress-induced separation Natural-
Stress-Moustakas2007 ; YAMANE20121 ; controlled spalling SPALLING-sapphire-
LEDbedell-sadana-2013 ; chemical etching of ZnO CLO-ZnO-Rogers2007 , Ga2O3
CLO-TsaiGa2o3-2011 , CrN CLO-CrN-2008-Yao , AlN CLO-AlN-2010-ChengFengLin  and
Nb2N nb2n-lift-off-2016  sacrificial layers; doping-selective electrochemical
etching CLO-nPlus-SiGaN-2009-Park ; substrate removal by grinding or etching
melnik_nikolaev_nikitina_vassilevski_dmitriev_1997 ; ETCHING-GaAs-
sumitomo-2001 ; Si-wafer-remove-Lesecq-2011 ; Mikulik-Si-etching-2013 ; void-
assisted separation USUI-VAS-2003 ; growth over patterned masks gogova2004elog
; HENNIG2008911-WSiN-ELOG ; lipski2010fabrication ; AMILUSIK201399 ; the use
of weakly bonded release layers like graphene Graphene-IBM-kim2014principle ,
BN BN-kobayashi2012layered , carbon CARBON-Phil-2016 ; and the use of
substrates with easy cleavage along the c-plane of a GaN epilayer Ga2O3-lift-
off-Gogova2012 ; SCAM-Matsuoka-2017 .

To separate a GaN film from a native GaN substrate several approaches were
proposed: chemical lift-off process InGaN-GaN-PEC-Youtsey2017 ; porous release
layers created by chemical LEE-KOH-POROUS-ETCH , electrochemical Mynbaeva-
porous-liftoff ; porous-lif-off-2013  or dry H2-porous-etching-YEH2011
etching of GaN substrate; controlled spalling SPALLING-bedell-sadana2017 ;
laser lift-off with an InGaN release layer LLO-InGaN-AKASAKI2014 ; ion
implantation Smart-Cut-GaN-Tauzin2005 ; Moutanabbir2010-ion-cleaving ; Ion-
Cut-Samsung-2013 ; huang2018defects ; and free-carrier-absorption laser lift-
off FCA-LLO-Virko2016 .

The lift-off methods based on the use of an intermediate layer require a
complicated epitaxy process. The ion implantation method does not allow
separating epitaxial device structures due to a large number of point defects
formed by the implantation process Ion-Cut-Samsung-2013 ; He-ion-damage-
HAN2017 ; huang2018defects . The free-carrier-absorption LLO method is limited
to separating undoped films from a heavily doped substrate FCA-LLO-Virko2016 .

We have proposed a laser slicing lift-off (LSLO) method for separating an
epitaxial structure from a bulk GaN substrate, one that does not require any
special release layers and can be used both for lifting off pure GaN films and
for lifting off GaN films with device structure shreter2013method . The method
is based on the effect of GaN decomposition induced by an ultra-short laser
pulse focused inside a bulk GaN material. The use of femtosecond lasers is
preferable as it allows to reduce the heat-affected zone and the corresponding
damage in the surrounding material stuart1996-ns-to-fs . This effect was
earlier used to create hollow 3D-channels inside bulk GaN FS-GaN-MACHINING-
MIDORIKAWA-2009  and to dice transparent wafers with laser stealth dicing
technology LSD-Kumagai-2007 .

The principle of the LSLO method is as follows. To lift-off a thin GaN film
with an InGaN multiple quantum well (MQW) device structure from a bulk GaN
substrate, the pulses of a near-infrared femtosecond laser are focused inside
the GaN layer, several microns under the surface (fig.1). The laser photon
energy is less than the InGaN bandgap, so both the InGaN device structure and
the GaN substrate are transparent at the wavelength of the laser. The laser
pulse energy is adjusted so that the nonlinear breakdown threshold is exceeded
only in the focus area and the GaN material in this area is decomposed.
Scanning of the focus position is performed in the $X$-$Y$-plane until a
continuous layer of decomposed material is created inside the wafer.

Figure 1: Experimental setup for the laser slicing lift-off process: a beam of
Yb-doped femtosecond laser with a pulse width of 350 fs and a pulse energy of
0.24 $\mu$J is focused inside a GaN wafer, 6 $\mu$m below the surface. The
decomposition of GaN takes place in the focus area of the laser beam. The
wafer is scanned in the $X$-$Y$ plane until a continuous layer of decomposed
material is created inside the wafer.

In this work, a proof-of-concept LED chip fabrication using the LSLO method is
described. First, a 5-$\mu$m thick GaN film with an InGaN LED structure was
grown homoepitaxially on a bulk GaN layer. After that, the film with the LED
structure was lifted-off from the substrate and transferred to a copper
carrier (fig. 2). The electroluminescence of the LED chip after the lift-off
was demonstrated. A detailed study of the LED chip for possible damage caused
by the LSLO process will be published separately.

A GaN-on-sapphire HVPE template was employed as a bulk-like substrate for the
InGaN LED structure epitaxy. The template was grown on a 2-inch c-plane (0001)
sapphire substrate using a two-stage growth process voronenkov2013two . The
thickness of the bulk-like HVPE GaN layer was 50 $\mu$m. The surface of the
HVPE GaN template was epi-polished by a chemical-mechanical polishing process.
After that, the InGaN MQW LED structure was deposited using a metal-organic
chemical vapor deposition system. The LED structure consisted of a 5-$\mu$m
thick $n^{+}$ doped GaN:Si ($n$=$4\cdot 10^{18}$ cm-3) layer; five pairs of
InGaN quantum wells with a well width of 2.3 nm and a barrier width of 11 nm;
and a 0.3-$\mu$m thick p-doped GaN:Mg layer (p=$3\cdot 10^{17}$ cm-3).

Figure 2: Schematic diagram of thin-film InGaN LED chip fabrication using the
laser slicing lift-off process.

The chip isolation and the LSLO process were performed with an ytterbium-
doped, diode-pumped solid-state laser with a wavelength of 1030 nm, a pulse
duration of 350 fs, and a pulse energy of 400 $\mu$J. A three-axis motorized
translation stage was used to perform scanning in the $X$-$Y$-plane and to
adjust the focus position along the $Z$ axis.

First, a 500$\mu$m$\times$500$\mu$m chip was isolated through the laser
scribing process – grooves with a depth of 6 $\mu$m and a width of 5 $\mu$m
were cut along the chip perimeter. Such grooves are necessary to release the
gas produced by the GaN decomposition during the LSLO process. A 20$\times$
dry objective was used to focus the laser beam on the surface of the GaN
layer. The pulse energy was 5 $\mu$J, pulse repetition rate was 1 kHz and the
scanning velocity was 1000 $\mu$m/s.

To perform the LSLO an oil immersion objective (100$\times$, NA = 1.25) with
an $n$=$1.518$ immersion oil was used. Before performing the LSLO, the
breakdown threshold was determined: the laser was focused at a depth of 6
$\mu$m below the surface and the laser pulse energy was gradually increased
until a breakdown in the focus area was observed at a pulse energy of 0.2
$\mu$J. After that, the LSLO was performed using a pulse energy of 0.24
$\mu$J. The laser beam, directed from the front side of the sample with a beam
focus located 6 $\mu$m under the surface was used to create a thin decomposed
layer 6 $\mu$m below the surface of the wafer. The scanning pattern exhibited
a set of parallel rows. The scanning velocity along the row was 1000 $\mu$m/s
and the pulse repetition rate was 2.5 kHz resulting in a pulse-to-pulse
distance along the row of 0.4 $\mu$m; the distance between the rows was 2
$\mu$m.

Figure 3: a) A photograph of the InGaN LED chip on a copper substrate with a
current probe connected, electroluminescence at 1 mA operating current. The
inset shows a microphotograph of this chip. b) SEM image of the Ga-face of the
chip after the LSLO process. c) SEM image of the N-face of the chip after the
LSLO process.

After the LSLO process, the LED chip was bonded to a temporary carrier wafer
using a cyanoacrylate adhesive. After the adhesive hardened, the chip was
debonded from the substrate. Then, the cyanoacrylate adhesive was dissolved,
and the chip was bonded to a copper substrate using an electrically conductive
silver epoxy as a bonding material.

To observe the electroluminescence of the chip, an electrical current probe
was connected directly to the surface of the chip, without applying any
current-spreading layers. The light emitting image of the LED chip on a copper
substrate is shown in fig. 3(a). The operating current was 1 mA. Blue-light
emission was observed under the probe, confirming that the MQW LED structure
was not destroyed during the LSLO process. A detailed study of the I-V
characteristics and the electroluminescence spectrum, as well as a detailed
comparison of the LED performance before and after the LSLO process will be
published separately.

The SEM image of the top side (Ga-face) of the LED chip after the LSLO process
is shown in fig. 3(b). No signs of damage caused by the LSLO were observed on
the surface. The SEM image of the bottom side (nitrogen-face) of the LED chip
after LSLO and debonding is shown in fig. 3(c). The remains of the decomposed
layer were observed. The surface roughness of the bottom side of the chip was
estimated using a tilted-view SEM image (not shown) to be about
1$\leavevmode\nobreak\ \mu$m, which corresponds to the focus depth of the
laser beam. The thickness of the separated film was measured using a side-view
SEM image (not shown) and was 5 $\mu$m, which corresponds to the position of
the focus below the surface.

Following is a brief discussion of the technological limits and possible
improvements of the method. The maximum thickness of the lifted-off layer is
limited by the defocusing of the laser beam caused by spherical aberration and
optical inhomogeneity of the material. Recent advances in the growth of pure
GaN SCIOCS-2017-pure  make it possible to obtain an optically homogeneous
material, and the adaptive optics techniques can compensate for spherical
aberration when focusing to a depth of hundreds of microns Salter-2014-depth .
The minimum thickness of the decomposed layer is limited by the depth of focus
of the laser beam, which can be reduced by the use of high-numerical-aperture
objectives, such as solid immersion lenses SIL-Kino-1999 . Processing
throughput is mainly limited by the scanning speed, therefore the use of
parallel multi-beam machining methods Kato-2005-microlens ;
Booth-2010-parallel  is highly desirable. Wafers with already deposited
metallization layers can also be processed. For this, the focused laser beam
is directed from the underside of the wafer. Other materials that decompose
irreversibly under the action of laser breakdown can also be processed. The
greatest practical interest of such materials is diamond, which is graphitized
by laser breakdown Konov-2008-diamond .

In conclusion, a 5-$\mu$m thick GaN film with the InGaN LED MQW structure was
lifted off from the bulk GaN layer by focusing near-infrared femtosecond laser
pulses 6 $\mu$m under the surface of the GaN wafer. The film was transferred
to a copper carrier and the electroluminescence of the LED structure was
demonstrated.

###### Acknowledgements.

The authors gratefully thank Dr. I. Guiney and Prof. C. Humphreys from
Cambridge University for the InGaN LED structure growth.

## References

  * (1) J. W. Chung, J. Lee, E. L. Piner, and T. Palacios, “Seamless On-Wafer Integration of Si(100) MOSFETs and GaN HEMTs,” _IEEE Electron Device Letters_ , vol. 30, no. 10, pp. 1015–1017, Oct 2009. [Online]. Available: https://doi.org/10.1109/LED.2009.2027914
  * (2) A. J. Trindade, B. Guilhabert, E. Y. Xie, R. Ferreira, J. J. D. McKendry, D. Zhu, N. Laurand, E. Gu, D. J. Wallis, I. M. Watson, C. J. Humphreys, and M. D. Dawson, “Heterogeneous integration of gallium nitride light-emitting diodes on diamond and silica by transfer printing,” _Opt. Express_ , vol. 23, no. 7, pp. 9329–9338, Apr 2015. [Online]. Available: http://doi.org/10.1364/OE.23.009329
  * (3) O. Moutanabbir and U. Gosele, “Heterogeneous Integration of Compound Semiconductors,” _Annual Review of Materials Research_ , vol. 40, no. 1, pp. 469–500, 2010. [Online]. Available: http://doi.org/10.1146/annurev-matsci-070909-104448
  * (4) J. W. Pomeroy, M. Bernardoni, D. C. Dumka, D. M. Fanning, and M. Kuball, “Low thermal resistance GaN-on-diamond transistors characterized by three-dimensional Raman thermography mapping,” _Applied Physics Letters_ , vol. 104, no. 8, p. 083513, 2014. [Online]. Available: https://doi.org/10.1063/1.4865583
  * (5) D. Iida, S. Kawai, N. Ema, T. Tsuchiya, M. Iwaya, T. Takeuchi, S. Kamiyama, and I. Akasaki, “Laser lift-off technique for freestanding GaN substrate using an In droplet formed by thermal decomposition of GaInN and its application to light-emitting diodes,” _Applied Physics Letters_ , vol. 105, no. 7, p. 072101, 2014. [Online]. Available: http://doi.org/10.1063/1.4893757
  * (6) M. K. Kelly, O. Ambacher, R. Dimitrov, R. Handschuh, and M. Stutzmann, “Optical Process for Liftoff of Group III-Nitride Films,” _physica status solidi (a)_ , vol. 159, no. 1, pp. R3–R4, 1996. [Online]. Available: http://doi.org/10.1002/1521-396X(199701)159:1%3CR3::AID-PSSA99993%3E3.0.CO;2-F
  * (7) A. D. Williams and T. Moustakas, “Formation of large-area freestanding gallium nitride substrates by natural stress-induced separation of GaN and sapphire,” _Journal of Crystal Growth_ , vol. 300, no. 1, pp. 37 – 41, 2007, first International Symposium on Growth of Nitrides. [Online]. Available: https://doi.org/10.1016/j.jcrysgro.2006.10.224
  * (8) K. Yamane, M. Ueno, H. Furuya, N. Okada, and K. Tadatomo, “Successful natural stress-induced separation of hydride vapor phase epitaxy-grown GaN layers on sapphire substrates,” _Journal of Crystal Growth_ , vol. 358, pp. 1 – 4, 2012. [Online]. Available: https://doi.org/10.1016/j.jcrysgro.2012.07.038
  * (9) S. W. Bedell, C. Bayram, K. Fogel, P. Lauro, J. Kiser, J. Ott, Y. Zhu, and D. Sadana, “Vertical Light-Emitting Diode Fabrication by Controlled Spalling,” _Applied Physics Express_ , vol. 6, no. 11, p. 112301, 2013. [Online]. Available: http://stacks.iop.org/1882-0786/6/i=11/a=112301
  * (10) D. J. Rogers, F. H. Teherani, A. Ougazzaden, S. Gautier, L. Divay, A. Lusson, O. Durand, F. Wyczisk, G. Garry, T. Monteiro, M. R. Correira, M. Peres, A. Neves, D. McGrouther, J. N. Chapman, and M. Razeghi, “Use of ZnO thin films as sacrificial templates for metal organic vapor phase epitaxy and chemical lift-off of GaN,” _Applied Physics Letters_ , vol. 91, no. 7, p. 071120, 2007. [Online]. Available: https://doi.org/10.1063/1.2770655
  * (11) T.-Y. Tsai, R.-H. Horng, D.-S. Wuu, S.-L. Ou, M.-T. Hung, and H.-H. Hsueh, “GaN Epilayer Grown on Ga2O3 Sacrificial Layer for Chemical Lift-Off Application,” _Electrochemical and Solid-State Letters_ , vol. 14, no. 11, pp. H434–H437, 2011. [Online]. Available: http://doi.org/10.1149/2.002111esl
  * (12) J. S. Ha, S. W. Lee, H. J. Lee, H. J. Lee, S. H. Lee, H. Goto, T. Kato, K. Fujii, M. W. Cho, and T. Yao, “The Fabrication of Vertical Light-Emitting Diodes Using Chemical Lift-Off Process,” _IEEE Photonics Technology Letters_ , vol. 20, no. 3, pp. 175–177, Feb 2008. [Online]. Available: http://doi.org/10.1109/LPT.2007.912491
  * (13) C.-F. Lin, J.-J. Dai, M.-S. Lin, K.-T. Chen, W.-C. Huang, C.-M. Lin, R.-H. Jiang, and Y.-C. Huang, “An AlN Sacrificial Buffer Layer Inserted into the GaN/Patterned Sapphire Substrate for a Chemical Lift-Off Process,” _Applied Physics Express_ , vol. 3, no. 3, p. 031001, 2010. [Online]. Available: https://doi.org/10.1143/APEX.3.031001
  * (14) D. J. Meyer, B. P. Downey, D. S. Katzer, N. Nepal, V. D. Wheeler, M. T. Hardy, T. J. Anderson, and D. F. Storm, “Epitaxial Lift-Off and Transfer of III-N Materials and Devices from SiC Substrates,” _IEEE Transactions on Semiconductor Manufacturing_ , vol. 29, no. 4, pp. 384–389, Nov 2016. [Online]. Available: http://doi.org/10.1109/TSM.2016.2599839
  * (15) J. Park, K. M. Song, S.-R. Jeon, J. H. Baek, and S.-W. Ryu, “Doping selective lateral electrochemical etching of GaN for chemical lift-off,” _Applied Physics Letters_ , vol. 94, no. 22, p. 221907, 2009. [Online]. Available: https://doi.org/10.1063/1.3153116
  * (16) Y. Melnik, A. Nikolaev, I. Nikitina, K. Vassilevski, and V. Dmitriev, “Properties of Free-Standing GaN Bulk Crystals Grown by HVPE,” _MRS Proceedings_ , vol. 482, p. 269, 1997. [Online]. Available: http://doi.org/10.1557/PROC-482-269
  * (17) K. Motoki, T. Okahisa, N. Matsumoto, M. Matsushima, H. Kimura, H. Kasai, K. Takemoto, K. Uematsu, T. Hirano, M. Nakayama, S. Nakahata, M. Ueno, D. Hara, Y. Kumagai, A. Koukitu, and H. Seki, “Preparation of Large Freestanding GaN Substrates by Hydride Vapor Phase Epitaxy Using GaAs as a Starting Substrate,” _Japanese Journal of Applied Physics_ , vol. 40, no. 2B, p. L140, 2001. [Online]. Available: https://doi.org/10.1143/JJAP.40.L140
  * (18) M. Lesecq, V. Hoel, A. L. des Etangs-Levallois, E. Pichonat, Y. Douvry, and J. C. D. Jaeger, “High Performance of AlGaN/GaN HEMTs Reported on Adhesive Flexible Tape,” _IEEE Electron Device Letters_ , vol. 32, no. 2, pp. 143–145, Feb 2011. [Online]. Available: http://doi.org/10.1109/LED.2010.2091251
  * (19) M. Lee, D. Mikulik, J. Kim, Y. Tak, J. Kim, M. Shim, Y. Park, U. Chung, E. Yoon, and S. Park, “A Novel Growth Method of Freestanding GaN Using In situ Removal of Si Substrate in Hydride Vapor Phase Epitaxy,” _Applied Physics Express_ , vol. 6, no. 12, p. 125502, 2013. [Online]. Available: https://doi.org/10.7567/APEX.6.125502
  * (20) Y. Oshima, T. Eri, M. Shibata, H. Sunakawa, K. Kobayashi, T. Ichihashi, and A. Usui, “Preparation of Freestanding GaN Wafers by Hydride Vapor Phase Epitaxy with Void-Assisted Separation,” _Japanese Journal of Applied Physics_ , vol. 42, no. 1A, p. L1, 2003. [Online]. Available: https://doi.org/10.1143/JJAP.42.L1
  * (21) D. Gogova, A. Kasic, H. Larsson, C. Hemmingsson, B. Monemar, F. Tuomisto, K. Saarinen, L. Dobos, B. Pécz, P. Gibart, and B. Beaumont, “Strain-free bulk-like GaN grown by hydride-vapor-phase-epitaxy on two-step epitaxial lateral overgrown GaN template,” _Journal of Applied Physics_ , vol. 96, no. 1, pp. 799–806, 2004. [Online]. Available: https://doi.org/10.1063/1.1753073
  * (22) C. Hennig, E. Richter, M. Weyers, and G. Trankle, “Freestanding 2-in GaN layers using lateral overgrowth with HVPE,” _Journal of Crystal Growth_ , vol. 310, no. 5, pp. 911 – 915, 2008, proceedings of the E-MRS Conference, Symposium G. [Online]. Available: https://doi.org/10.1016/j.jcrysgro.2007.11.102
  * (23) F. Lipski, T. Wunderer, S. Schwaiger, and F. Scholz, “Fabrication of freestanding 2”-GaN wafers by hydride vapour phase epitaxy and self-separation during cooldown,” _Physica Status Solidi (a)_ , vol. 207, no. 6, pp. 1287–1291, 2010. [Online]. Available: https://doi.org/10.1002/pssa.200983517
  * (24) M. Amilusik, T. Sochacki, B. Lucznik, M. Bockowski, B. Sadovyi, A. Presz, I. Dziecielewski, and I. Grzegory, “Analysis of self-lift-off process during HVPE growth of GaN on MOCVD-GaN/sapphire substrates with photolitographically patterned Ti mask,” _Journal of Crystal Growth_ , vol. 380, pp. 99 – 105, 2013. [Online]. Available: https://doi.org/10.1016/j.jcrysgro.2013.06.005
  * (25) J. Kim, C. Bayram, H. Park, C.-W. Cheng, C. Dimitrakopoulos, J. A. Ott, K. B. Reuter, S. W. Bedell, and D. K. Sadana, “Principle of direct van der Waals epitaxy of single-crystalline films on epitaxial graphene,” _Nature communications_ , vol. 5, p. 4836, 2014. [Online]. Available: https://doi.org/10.1038/ncomms5836
  * (26) Y. Kobayashi, K. Kumakura, T. Akasaka, and T. Makimoto, “Layered boron nitride as a release layer for mechanical transfer of GaN-based devices,” _Nature_ , vol. 484, no. 7393, p. 223, 2012. [Online]. Available: https://doi.org/10.1038/nature10970
  * (27) A. S. Altakhov, R. I. Gorbunov, L. A. Kasharina, F. E. Latyshev, V. A. Tarala, and Y. G. Shreter, “Amorphous carbon buffer layers for separating free gallium nitride films,” _Technical Physics Letters_ , vol. 42, no. 11, pp. 1076–1078, Nov 2016. [Online]. Available: https://doi.org/10.1134/S106378501611002X
  * (28) K. Kachel, M. Korytov, D. Gogova, Z. Galazka, M. Albrecht, R. Zwierz, D. Siche, S. Golka, A. Kwasniewski, M. Schmidbauer, and R. Fornari, “A new approach to free-standing GaN using $\beta$-Ga2O3 as a substrate,” _CrystEngComm_ , vol. 14, pp. 8536–8540, 2012. [Online]. Available: http://dx.doi.org/10.1039/C2CE25976A
  * (29) K. Ohnishi, M. Kanoh, T. Tanikawa, S. Kuboya, T. Mukai, and T. Matsuoka, “Halide vapor phase epitaxy of thick GaN films on ScAlMgO 4 substrates and their self-separation for fabricating freestanding wafers,” _Applied Physics Express_ , vol. 10, no. 10, p. 101001, 2017. [Online]. Available: https://doi.org/10.7567/APEX.10.101001
  * (30) Y. Chris, M. Robert, R. Rekha, F. Kamran, X. Andy, B. Ed, W. Jingshan, F. Patrick, C. Timothy, C. Eric, and G. Louis, “Wafer-scale epitaxial lift-off of GaN using bandgap-selective photoenhanced wet etching,” _physica status solidi (b)_ , vol. 254, no. 8, p. 1600774, 2017. [Online]. Available: https://doi.org/10.1002/pssb.201600774
  * (31) C.-Y. Lee, Y.-P. Lan, P.-M. Tu, S.-C. Hsu, C.-C. Lin, H.-C. Kuo, G.-C. Chi, and C.-Y. Chang, “Natural substrate lift-off technique for vertical light-emitting diodes,” _Applied Physics Express_ , vol. 7, no. 4, p. 042103, 2014. [Online]. Available: https://doi.org/10.7567/APEX.7.042103
  * (32) M. G. Mynbaeva, A. E. Nikolaev, A. A. Sitnikova, and K. D. Mynbaev, “HVPE homo-epitaxial growth of GaN on porous substrates,” _CrystEngComm_ , vol. 15, pp. 3640–3646, 2013. [Online]. Available: http://dx.doi.org/10.1039/C3CE27099H
  * (33) S. Huang, Y. Zhang, B. Leung, G. Yuan, G. Wang, H. Jiang, Y. Fan, Q. Sun, J. Wang, K. Xu, and J. Han, “Mechanical Properties of Nanoporous GaN and Its Application for Separation and Transfer of GaN Thin Films,” _ACS Applied Materials and Interfaces_ , vol. 5, no. 21, pp. 11 074–11 079, 2013, pMID: 24125198. [Online]. Available: http://doi.org/10.1021/am4032345
  * (34) Y.-H. Yeh, K.-M. Chen, Y.-H. Wu, Y.-C. Hsu, T.-Y. Yu, and W.-I. Lee, “Hydrogen etching of GaN and its application to produce free-standing GaN thick films,” _Journal of Crystal Growth_ , vol. 333, no. 1, pp. 16 – 19, 2011. [Online]. Available: https://doi.org/10.1016/j.jcrysgro.2011.08.022
  * (35) S. W. Bedell, P. Lauro, J. A. Ott, K. Fogel, and D. K. Sadana, “Layer transfer of bulk gallium nitride by controlled spalling,” _Journal of Applied Physics_ , vol. 122, no. 2, p. 025103, 2017. [Online]. Available: https://doi.org/10.1063/1.4986646
  * (36) A. Tauzin, T. Akatsu, M. Rabarot, J. Dechamp, M. Zussy, H. Moriceau, J. F. Michaud, A. M. Charvet, L. D. Cioccio, F. Fournel, J. Garrione, B. Faure, F. Letertre, and N. Kernevez, “Transfers of 2-inch GaN films onto sapphire substrates using Smart CutTM technology,” _Electronics Letters_ , vol. 41, no. 11, pp. 668–670, May 2005. [Online]. Available: http://doi.org/10.1049/el:20051038
  * (37) O. Moutanabbir and U. Gösele, “Bulk GaN Ion Cleaving,” _Journal of Electronic Materials_ , vol. 39, no. 5, pp. 482–488, May 2010. [Online]. Available: https://doi.org/10.1007/s11664-010-1100-0
  * (38) R. B.-K. Chung, D. Kim, S.-K. Lim, J.-S. Choi, K.-J. Kim, B.-H. Lee, K. S. Jung, H.-J. Kim-Lee, W. J. Lee, B. Park, and K. Woo, “Layer-Transferred GaN Template by Ion Cut for Nitride-Based Light-Emitting Diodes,” _Applied Physics Express_ , vol. 6, no. 11, p. 111005, 2013. [Online]. Available: https://doi.org/10.7567/APEX.6.111005
  * (39) K. Huang, T. You, Q. Jia, A. Yi, S. Zhang, R. Zhang, J. Lin, M. Zhou, W. Yu, B. Zhang _et al._ , “Defects induced by MeV H+ implantation for exfoliating of free-standing GaN film,” _Applied Physics A_ , vol. 124, no. 2, p. 118, 2018. [Online]. Available: https://doi.org/10.1007/s00339-017-1508-y
  * (40) M. V. Virko, V. S. Kogotkov, A. A. Leonidov, V. V. Voronenkov, Y. T. Rebane, A. S. Zubrilov, R. I. Gorbunov, P. E. Latyshev, N. I. Bochkareva, Y. S. Lelikov, D. V. Tarhin, A. N. Smirnov, V. Y. Davydov, and Y. G. Shreter, “On the laser detachment of n-GaN films from substrates, based on the strong absorption of IR light by free charge carriers in n+-GaN substrates,” _Semiconductors_ , vol. 50, no. 5, pp. 699–704, May 2016. [Online]. Available: https://doi.org/10.1134/S1063782616050250
  * (41) Y. Han, J. Peng, B. Li, Z. Wang, K. Wei, T. Shen, J. Sun, L. Zhang, C. Yao, N. Gao, X. Gao, L. Pang, Y. Zhu, H. Chang, M. Cui, P. Luo, Y. Sheng, H. Zhang, L. Zhang, X. Fang, S. Zhao, J. Jin, Y. Huang, C. Liu, P. Tai, D. Wang, and W. He, “Lattice disorder produced in GaN by He-ion implantation,” _Nuclear Instruments and Methods in Physics Research Section B: Beam Interactions with Materials and Atoms_ , vol. 406, pp. 543 – 547, 2017, the 27th International Conference on Atomic Collisions in Solids. [Online]. Available: https://doi.org/10.1016/j.nimb.2016.12.039
  * (42) Y. G. Shreter, Y. T. Rebane, and A. V. Mironov, “Method of Separating Surface Layer of Semiconductor Crystal Using a Laser Beam Perpendicular to the Separating Plane,” Sep. 26 2013, US Patent App. 13/990,359.
  * (43) B. C. Stuart, M. D. Feit, S. Herman, A. M. Rubenchik, B. W. Shore, and M. D. Perry, “Nanosecond-to-femtosecond laser-induced breakdown in dielectrics,” _Phys. Rev. B_ , vol. 53, pp. 1749–1761, Jan 1996. [Online]. Available: https://doi.org/10.1103/PhysRevB.53.1749
  * (44) S. Nakashima, K. Sugioka, and K. Midorikawa, “Fabrication of microchannels in single-crystal GaN by wet-chemical-assisted femtosecond-laser ablation,” _Applied Surface Science_ , vol. 255, no. 24, pp. 9770 – 9774, 2009, proceedings of the Sixth International Conference on Photo-Excited Processes and Applications (6-ICPEPA). [Online]. Available: https://doi.org/10.1016/j.apsusc.2009.04.159
  * (45) M. Kumagai, N. Uchiyama, E. Ohmura, R. Sugiura, K. Atsumi, and K. Fukumitsu, “Advanced Dicing Technology for Semiconductor Wafer Stealth Dicing,” _IEEE Transactions on Semiconductor Manufacturing_ , vol. 20, no. 3, pp. 259–265, Aug 2007. [Online]. Available: http://doi.org/10.1109/TSM.2007.901849
  * (46) V. V. Voronenkov, N. I. Bochkareva, R. I. Gorbunov, P. E. Latyshev, Y. S. Lelikov, Y. T. Rebane, A. I. Tsyuk, A. S. Zubrilov, U. W. Popp, M. Strafela, and Y. Shreter, “Two modes of HVPE growth of GaN and related macrodefects,” _Physica Status Solidi (c)_ , vol. 10, no. 3, pp. 468–471, 2013. [Online]. Available: https://doi.org/10.1002/pssc.201200701
  * (47) H. Fujikura, T. Konno, T. Yoshida, and F. Horikiri, “Hydride-vapor-phase epitaxial growth of highly pure GaN layers with smooth as-grown surfaces on freestanding GaN substrates,” _Japanese Journal of Applied Physics_ , vol. 56, no. 8, p. 085503, 2017. [Online]. Available: https://doi.org/10.7567/JJAP.56.085503
  * (48) P. S. Salter, M. Baum, I. Alexeev, M. Schmidt, and M. J. Booth, “Exploring the depth range for three-dimensional laser machining with aberration correction,” _Opt. Express_ , vol. 22, no. 15, pp. 17 644–17 656, Jul 2014. [Online]. Available: http://doi.org/10.1364/OE.22.017644
  * (49) G. S. Kino, “Applications and theory of the solid immersion lens,” _Proc.SPIE_ , vol. 3609, pp. 3609 – 3609 – 11, 1999. [Online]. Available: https://doi.org/10.1117/12.351043
  * (50) J.-i. Kato, N. Takeyasu, Y. Adachi, H.-B. Sun, and S. Kawata, “Multiple-spot parallel processing for laser micronanofabrication,” _Applied Physics Letters_ , vol. 86, no. 4, p. 044102, 2005. [Online]. Available: https://doi.org/10.1063/1.1855404
  * (51) A. Jesacher and M. J. Booth, “Parallel direct laser writing in three dimensions with spatially dependent aberration correction,” _Opt. Express_ , vol. 18, no. 20, pp. 21 090–21 099, Sep 2010. [Online]. Available: http://doi.org/10.1364/OE.18.021090
  * (52) T. Kononenko, M. Meier, M. Komlenok, S. Pimenov, V. Romano, V. Pashinin, and V. Konov, “Microstructuring of diamond bulk by IR femtosecond laser pulses,” _Applied Physics A_ , vol. 90, no. 4, pp. 645–651, Mar 2008. [Online]. Available: https://doi.org/10.1007/s00339-007-4350-9


∎

11institutetext: Barcelona Supercomputing Center (BSC)
11email<EMAIL_ADDRESS>
_This work has been submitted to the Springer for possible publication.
Copyright may be transferred without notice, after which this version may no
longer be accessible._

# cuConv: A CUDA Implementation of Convolution for CNN Inference

Marc Jorda  Pedro Valero-Lara  Antonio J. Peña

(Received: date / Accepted: date)

###### Abstract

Convolutions are the core operation of deep learning applications based on
Convolutional Neural Networks (CNNs). Current GPU architectures are highly
efficient for training and deploying deep CNNs, and hence, these are largely
used in production for this purpose. State–of–the–art implementations,
however, present a lack of efficiency for some commonly used network
configurations.

In this paper we propose a GPU-based implementation of the convolution
operation for CNN inference that favors coalesced accesses, without requiring
prior data transformations. Our experiments demonstrate that our proposal
yields notable performance improvements in a range of common CNN forward
propagation convolution configurations, with speedups of up to 2.29$\times$
with respect to the best implementation of convolution in cuDNN, hence
covering a relevant region in currently existing approaches.

###### Keywords:

Coalescing Convolutional neural networks cuDNN Deep learning GPU convolution

##  1 Introduction

Even though Convolutional Neural Networks (CNNs)cognitron ; neocognitron ;
lecun98  were presented many years ago, since the appearance of CNNs in the
ImageNet Large Scale Visual Recognition Challenge (ILSVRC)
NIPS2012_4824_alexnet ; ILSVRC15  in 2012, CNNs and other deep neural network
architectures have received considerable more attention than in previous
years. This interest comes from their impressive results in tasks like image
classification, speech recognition, and natural language processing
GuWKMSSLWW15_cnns . Two main factors were necessary to enable the success of
CNNs: the availability of very large datasets of manually annotated input
data, and the high performance of current computing systems. Large datasets
are needed to train the deep neural network parameters until a highly accurate
result is reached. In turn, such amount of data required the use of high
performance computing accelerators to keep the training time of deep neural
networks within reasonable limits. GPUs are probably the most common type of
accelerator used for deep learning, among others like FPGAs or custom
solutions, and all of the main deep learning frameworks currently support
execution on GPUs.

High performance computing is also important when CNN-based applications are
deployed. Short response times are one of the most relevant parameters in
terms of user satisfaction. Moreover, short latency requirements are mandatory
for applications where delays in the response time pose safety implications.
Deep learning-based safety-critical applications, like autonomous vehicles,
will attain widespread use in the near future, which highlights the need for
high performance CNN inference.

CNNs differ from other deep neural networks due to the inclusion of
_convolutional layers_. The outputs of these layers are the result of weighted
sums of inputs, like in fully-connected layers. The difference strives in the
inputs that are involved in the computation of each output, and the fact that
convolutional layers feature a smaller number of weights which are reused for
several outputs. Outputs of a convolutional layer are not a function of all
the inputs as in fully-connected layers; they only depend on a contiguous
subset of the inputs (known as receptive field). This considerably reduces the
computation cost of these layers, when compared to fully-connected layers. It
is also interesting that, in convolutional layers, computational and storage
costs no longer depend on the number of inputs and outputs. The computational
cost depends on the size of the contiguous subset of inputs, while the reuse
of weights reduces the memory required to store them. To implement the
weighted sums, convolutional layers use an operation called _convolution_ ,
which gives them the name.

Most of the execution time of a convolutional layer is spent performing
convolutions. Since most of the layers in CNNs are convolutional, convolutions
account for a large part of the overall network execution time. Several works
target the optimization of convolution operations for GPUs. Some perform data
transformations to be able to exploit already existing high-performance
functions like GEMM Chetlur_cuDNN14 . Other works rely on algorithmic
optimizations to reduce the computational cost of convolutions, thus improving
their performance Lavin15b_winograd ; VasilacheJMCPL14_fft .

In this paper, we present the design and implementation of a GPU convolution
algorithm for NVIDIA GPUs. Our approach is based on efficiently exploiting the
GPU execution resources and in-core memories. Our design also does an
efficient use of the GPU memory bandwidth, performing coalesced accesses
without the need for costly data transformations before the main computation
kernel. To evaluate our implementation of convolution for GPUs, we compare its
performance with all the GPU convolution algorithms provided by the cuDNN
library. This library is developed by NVIDIA and contains several
implementations of convolution based on the current state–of–the–art
algorithms. Currently, the convolutions and other deep learning operations
provided by cuDNN are used as the GPU backend in most of the widely used deep
learning frameworks, like PyTorch pytorch  or TensorFlow tensorflow . For the
experimental evaluation, we use the forward propagation convolution parameter
configurations from five of the most well known CNN architectures, namely
AlexNet NIPS2012_4824_alexnet , GoogleNet SzegedyLJSRAEVR14_googlenet ,
ResNet-50 HeZRS15_resnet , SqueezeNet IandolaMAHDK16_squeezenet , and VGG19
SimonyanZ14a_vgg . This ensures that our evaluation only focuses on relevant
convolution configurations that are actually present in CNNs used in deep
learning research and CNN-based applications. In total, we test more than 600
different convolution parameter configurations. As part of the evaluation, we
also profile the execution of some representative parameter configurations to
further explain the observed performance results.

This paper is organized as follows. In Section 2 we introduce the concepts
required to understand this work. In Section 3 we present the design of our
GPU convolution implementation and discuss some of the trade-offs we faced
during the design. In Section 4 we describe the experimental setup and present
the obtained performance results. In Section 5 we discuss other work related
to convolution optimizations. In Section 6 we summarize this work’s main
conclusions and present possible future directions to improve our design.

##  2 Background

In this section we provide background information that will help readers
understand the rest of the paper. First we present some details on
convolutional layers and the convolution operation itself. Next, we introduce
a few GPU architecture and programming model terms and naming conventions that
will be used in the following sections. Finally, we briefly describe the three
major state–of–the–art convolution algorithms, highlighting their strengths
and weaknesses.

###  2.1 Convolutional Layers

Convolutional layers in CNNs receive a set of _N_ 3D inputs (the _input
batch_) and generate an equally sized set of 3D outputs. The elements of each
output volume are computed as a weighted sum of some of the elements of the
corresponding input, which is usually transformed afterwards by a non-linear
function (like sigmoid or ReLU abs-1710-05941_activations ). These layers use
convolution operations to perform the weighted sum of the input elements,
hence their name. The weights of the layer are represented as a set of _M_
volumes called _filters_ or _kernels_. Since the term _kernel_ is also used to
describe a GPU function, we will use _filters_ to refer to the weights of
convolutional layers. The convolution of a 3D input with a given filter
produces a matrix. The matrices from the convolutions of a given input with
all the layer’s filters are stacked in the Z dimension to produce the
corresponding output (see Figure 1). This process is repeated for all inputs
to obtain the final set of _N_ output volumes.

Figure 1: Convolution operations in a convolutional layer. The matrix produced
by the convolution of Input 0 with Filter 0 is highlighted in light blue. The
darker output element is the result of the dot product of Filter 0 with the
highlighted subvolume of Input 0.

A _convolution operation_ in the context of convolutional layers is a 2D
discrete convolution, with the particularity that it uses a 3D input and
filter. Each output element is computed as the dot product of the filter with
a subvolume of the input, as depicted in Figure 1. The input subvolume used
for each output element is different and depends on its coordinates. The top
left element uses the subvolume in the top left corner of the input, and the
rest of output elements are computed by translating the subvolume over the X
and Y dimensions of the input (there is no translation over the Z dimension).
The size of the input along the Z dimension (referred as the _number of
channels_) and the size of the filter along the Z dimension must be equal, for
the dot products to be well defined. In this paper, we refer to the size in
the Z dimension of inputs/filters/outputs as their _depth_. The translation
distance among two contiguous positions is called _stride_ , and defines the
width and height of the output. Increasing the stride reduces the number of
positions the filter can translate over the input, resulting in smaller
outputs. _Padding_ is another parameter related to the output width and
height. It defines the number of padding rows and columns added to the input
along the X and Y dimensions, to increase its size. This increases the
positions for the filter translation, increasing the output width and height.
In CNNs, it is common to use a padding of $\frac{W_{f}-1}{2}$ elements in each
side along the X dimension, and a padding of $\frac{H_{f}-1}{2}$ elements
along the Y dimension. Combined with a stride of 1 element in both dimensions,
these values generate an output with the same width and height of the input.

Inputs, filters and outputs are all sets of 3D volumes. In CNNs, these sets
are represented as 4D arrays usually called _tensors_. In terms of the memory
layout used to store the tensors, there is no standard agreement across the
different deep learning frameworks. Several orderings of the four dimensions
are being used, but two of the most common layouts are NCHW and CHWN
mem_layout_Li16 . These letters stand for height (H), width (W), channels (C,
another name for _depth_), and N as the dimension representing the index of
the volume within the set. The fourth dimension in the abbreviations is that
with the elements contiguous in memory.

Most of the current deep learning frameworks use the operations from the cuDNN
library as their GPU implementation of convolutions and other deep learning
operations. cuDNN provides several implementations of the convolution
operation, with different performance depending on the convolution parameters
(input and filter sizes, stride, batch size, etc.). Because of this, several
frameworks perform an initial exploration to choose the best-performing
implementation of convolution for each convolutional layer.

###  2.2 GPU Architecture and Naming Conventions

Figure 2: Schematic of our target GPU architecture.

Throughout this paper we use NVIDIA terminology to refer to GPU architecture
elements and to concepts from the GPU programming and execution models. A
kernel is a function that contains the code that will be executed by each GPU
thread. The number of GPU threads that will execute a kernel is specified when
the kernel is launched, and is referred as the thread grid (from 1 to 3
dimensions). Threads are grouped into thread blocks. Threads in the same
thread block are scheduled together in the same Streaming Multiprocessor (SM),
and may be synchronized through thread barriers. Thread blocks are independent
of each other and cannot synchronize their execution. The GPU execution model
does not provide any guarantee on the thread blocks execution order; all
possible interleavings of thread block execution are valid. This provides
adaptive parallelism: a kernel launch with enough thread blocks may exploit
the parallelism provided by GPUs with different number of SMs without
requiring any modification from the developer.

In terms of the memory subsystem of our target GPU architecture, depicted in
Figure 2, we assume a non-coherent 2-level cache hierarchy with a private L1
cache per SM and an L2 cache shared by all the SMs. This configuration is
commonly found in discrete NVIDIA GPUs to connect the execution cores to the
on-board HBM2 memory. Besides the L1 cache, there are two other types of in-
core memories. Each SM features a software-managed scratchpad memory, called
_shared memory_ in CUDA terminology. Each thread block running in an SM
features its own space within the shared memory. Only threads of the same
thread block may share data through shared memory. Accesses to shared memory
are one or two orders of magnitude faster than accesses to global memory (that
miss/bypass the L1 cache) because these do not have to travel across the
shared on-chip interconnect to reach the L2, nor reach global memory in case
of an L2 miss. The last in-core memory is the register file, which stores the
kernel’s local variables, which are private to each thread.

Exploiting in-core memories is one of the most important aspects to attain a
high-performance GPU code, similar to a proper workload management to maximize
cache hits in the CPU case. However, for GPUs, it is also important to
consider the memory access pattern of contiguous threads. Execution units in
GPU cores are vector units of a certain width, and a _warp_ describes the
group of GPU threads that execute together as a single vector instruction. For
memory loads and stores, threads in a warp may access any arbitrary address,
but since these are executed in lock-step, all the threads will be stalled
until all the memory transactions are fulfilled. Thus, to attain the minimum
memory instruction cost, all threads in a warp have to access addresses within
the same cache line (usually referred as _coalesced_ accesses). Current NVIDIA
GPUs have a warp width of 32 threads, and an L1 cache line size of 128 bytes.

###  2.3 Convolution Algorithms

There are various algorithms that may be used to compute a convolution. The
first option is to directly apply the convolution formula, but there are more
complex approaches. The convolution inputs and filters can be transformed to
be able to exploit high-performance matrix-matrix multiply (GEMM)
implementations Chetlur_cuDNN14 . Winograd’s minimal filtering algorithm
Lavin15b_winograd  and Fast Fourier Transform (FFT) MathieuHL13_fft ;
VasilacheJMCPL14_fft  are two other approaches to compute convolutions, both
based on reducing the number of multiplications (but increasing the number of
additions) to reduce the computational cost of the computation.

####  2.3.1 GEMM

The GEMM-based algorithm transforms the filters into an intermediate matrix,
where each row contains the elements of a given filter, flattened as a 1D
array. The inputs are also transformed to an intermediate matrix, in such a
way that the dot products computed by the GEMM of these two matrices match the
ones that would be computed by the convolution formula. Finally, the GEMM’s
output matrix is transformed to obtain the convolution output with the
appropriate format. GEMM-based algorithms may benefit from already existing
high-performance GEMM on GPUs DongarraICCS17 . However, the intermediate
matrices require a significant amount of memory. Specifically, the one derived
from the convolution inputs has to store duplicate elements due to the overlap
of the filter positions when the stride is smaller than the filter, which
require more memory than the inputs themselves. To avoid this memory cost,
some variants of the GEMM based convolution algorithms do not explicitly store
the intermediate matrices and do the transformation when loading the data in
the in-core memories (e.g. shared memory, registers).

####  2.3.2 Winograd

The Winograd-based algorithm applies the transformations described in
winograd1980arithmetic  to minimize the number of multiplications needed for
the computation of the convolution. This number is reduced by increasing the
number of additions, which are usually less costly to compute. The number of
additions increases quadratically with the size of the input, making this
technique better suited for small inputs. For large inputs, the Winograd-based
algorithm divides the input in tiles of $p\times{}p$ elements, which have to
overlap by $p-2$ elements. After computing the minimal convolution for each
tile, the intermediate results are combined to obtain the final convolution
output.

####  2.3.3 FFT

The FFT-based convolution algorithms exploit the property that the convolution
in the time domain is equal to point-wise multiplication in the Fourier
(frequency) domain. The algorithm computes the FFT of the convolution inputs,
then performs the point-wise multiplication followed by an inverse FFT to get
the convolution output. While these transformations may be too costly to
compute a single convolution, their cost can be amortized by several
convolutions. This is the case of convolutional layers where the FFT of each
input and filter can be reused for all the convolutions where these
participate. Thus, the potential improvement of FFT-based algorithms increases
with larger number of inputs and/or larger number of filters.

##  3 Design and Implementation

In this section we present the details of the design of our GPU-based
algorithm to compute convolution. We also discuss and pay particular attention
to some relevant implementation trade-offs.

As introduced earlier, our efforts focus on efficiently exploiting the GPU
hardware features to implement a fast GPU implementation for convolutional
operations. To implement a fast GPU kernel, an effective exploitation of the
top-level hierarchies of memory (shared memory and registers) is of vital
importance, in order to maximize the reuse of data. This is in fact
particularly necessary for those throughput-oriented architectures, such as
GPUs, where the accesses to memory are impacted by a high latency. In the
computation performed by a convolutional layer we find two levels of data
reuse:

  1. 1.

At the layer level. Convolutional layers receive a set of inputs that have to
be convolved with a set of filters. At this level we may exploit the reuse of
data in both, the filters and the inputs. Each input is convolved with all the
filters or, the other way around, each filter is convolved with all the
inputs.

  2. 2.

At the convolution level. Inside a given convolution for a specific pair of
input and filter, the filter elements (i.e. weights) are shared by several
elements of the same input, which is one of the main differences when compared
to the fully-connected layers. As the filter is passed (computed) through all
the elements of the input, each filter row along the Z dimension is point-wise
multiplied with the corresponding rows along the Z dimension of the input (to
make the presentation less verbose, for the rest of this section we will refer
to the rows along the Z axis simply as _rows_). Besides, most of the
convolutional operations to be performed present a stride of 1 in both
dimensions, X (horizontal) and Y (vertical), which helps to a higher data
reuse, since the input elements which share one row of a given filter generate
a cube composed by contiguous rows (memory positions), as depicted in Figure
3. This scenario is particularly well-suited for GPU computing.

Figure 3: Reuse of input data for two example rows of a filter (highlighted in
blue and orange), for a convolution with a stride of 1. The highlighted input
rows (two sets of 6x6 rows) are point-wise multiplied with the filter row of
the same color during the convolution computation.

Since the filter size is usually smaller than the input size, it is better to
exploit the reuse of data from the filter side than from the input side,
because it requires less shared memory and each value brought to the shared
memory will be reused more times. It is important to note that the reuse of
data for those rows closer to the edges involves less filter-input operations.
This is particularly negative in the rows positioned in the corners, which may
be only used once. On the other extreme, those rows positioned around the
center of the input (with coordinates $W_{f}-1\leq x\leq W-W_{f}$ and
$H_{f}-1\leq y\leq H-H_{f}$) present the best ratio of data reuse ($H_{f}\cdot
W_{f}$ times per element). The data reuse of any row of a particular filter is
equal to $(W-W_{f}+1)\cdot(H-H_{f}+1)$, which is the product of the height and
width of the outputs of the convolution.

Our GPU convolution is divided into two stages, implemented in two separate
kernels. In the first stage, we compute all the dot products between input and
filter rows that are needed by all the convolutions of the convolutional
layer. The result of this first step is a set of $H_{f}\cdot W_{f}\cdot N\cdot
M$ temporary matrices, each of size $H-H_{f}+1$ by $W-W_{f}+1$, as depicted by
the scheme from Figure 4. In other words, the first kernel applies a reduction
to the contiguous input 3D patches along the Z dimension for all filters and
inputs, resulting in a collection of temporary 2D matrices with the size of an
output X-Y plane.

The second stage collects and sums the $H_{f}\cdot W_{f}$ temporary matrices
resulting from a given input-filter pair, and stores the result as the
corresponding output plane. In Figure 4, the output plane for the second stage
is the result of adding the matrices of partial results depicted in the
scheme. Output elements are the sum one element from each 2D matrix (i.e.
filter row), which is one of the reasons why this final reduction is performed
in a separate kernel.

The division into two stages enables us to exploit the reuse of filter
elements during the first stage. Each thread block is responsible for the dot
products of one (possibly more) filter row with all the rows from the inputs
that interact with it in any of the convolutions. In order to leverage the in-
core memory regions, the threads of a given thread block collaborate to bring
the filter row into a shared memory buffer. Once the filter data is ready,
each thread computes one (or several) of the dot products involving the
assigned filter row, and stores the result in a temporary output array.

Depending on the input size and the number of inputs and filters, the dot
products of one filter row may be split into several thread blocks. This may
be necessary when the number of threads per block required to perform all the
dot products is larger than the maximum number of threads per block allowed by
the GPU. It can also be a possible optimization to increase the thread block-
level parallelism and exploit all the SMs in configurations that would
otherwise feature a low number of thread blocks. However, this increases the
overall amount of long-latency memory accesses because the kernel has to move
the same filter row data several times from global memory to shared memory,
which may overcome the possible benefit of increasing the thread block-level
parallelism (the GPU L2 cache may reduce this impact depending on workload
size and the temporal locality of the thread blocks which access the same
filter row).

Figure 4: Stages of our implementation of convolution depicted for an
arbitrary input and the first filter of the convolutional layer. Scalar
products in the first stage generate partial results which are aggregated in
the second stage to obtain the final output elements.

In terms of data layout, our implementation supports tensors in NCHW format,
which is one of the most common memory layouts in CNNs mem_layout_Li16 . In
this format, elements in the X dimension are stored in contiguous memory
locations. Then it stores the rows along Y, the planes along Z, and finally
the rest of the 3D inputs/filter/outputs as the fourth dimension. The
advantage of this format, given the fact that the input rows that reuse a
given filter row in a stride 1 convolution are contiguous in the X dimension,
is that the reads of the input data performed by the first-stage kernel are
almost perfectly coalesced. Unlike other implementations, like the GEMM-based
convolution Chetlur_cuDNN14 , the coalesced accesses are attained without the
need for data-layout transformations of the inputs to adapt their layout to
the memory access pattern of the kernel. For convolution configurations where
the width of the contiguous input 3D patch read by a thread block is multiple
of 32, all the warps will perform coalesced accesses. For other
configurations, some of the warps will read elements from two different rows
along the X axis of the 3D patch, which introduces a gap of $W_{f}-1$ elements
in between the memory positions accessed by two contiguous threads (see Figure
3). In such a case, the memory load instruction will issue at least 2 memory
accesses, instead of the ideal single access. However, since this only occurs
to some of the warps, the cost of transforming the input memory layout to
always feature perfectly coalesced reads would overcome the benefits obtained
by avoiding the additional memory accesses.

For convolutions which involve filters of size 1$\times$1, the second kernel
is not necessary because each output element only depends on one input and one
filter row, and thus the outputs of the first kernel are already the final
output elements. In this case, the first-stage kernel stores its outputs with
the proper memory layout to represent the final 4D output tensor, avoiding the
need to execute the second-stage kernel.

##  4 Experimental Evaluation

Following the study of access , to analyze the performance of our proposed
implementation of convolution we selected all the forward propagation
convolutional layer configurations from five widely known CNNs: AlexNet
NIPS2012_4824_alexnet , GoogleNet SzegedyLJSRAEVR14_googlenet , ResNet-50
HeZRS15_resnet , SqueezeNet IandolaMAHDK16_squeezenet , and VGG19
SimonyanZ14a_vgg , identifying the input and filter sizes, the number of
filters, and the depth (see Table 1). For the number of inputs (batch size),
we evaluate a range of characteristic values (1, 8, 16, 32, 64, 128, 256)
because this parameter is usually selected by the end user of the CNN. The
final set includes configurations with the following filter sizes: 1$\times$1,
3$\times$3, and 5$\times$5\. The most common filter size is 1$\times$1, which
features a particular purpose within CNNs. Convolutional layers with this
filter size are used as a means to control the computational cost of the
following convolutional layer, which usually has a larger filter size. Tuning
the number of filters in the 1$\times$1 layer, to reduce the depth of the
outputs, directly translates in less computational load in the following
layer.

Table 1: Summary of features of the convolution operations (with stride 1) present in the five well-known CNNs we selected for our evaluation access . | GoogleNet | SqueezeNet | AlexNet | Resnet50 | VGG19
---|---|---|---|---|---
# Distinct convolution configurations | 42 | 21 | 4 | 12 | 9
Network input size |  224$\times$224$\times$3 |  224$\times$224$\times$3 |  224$\times$224$\times$3 |  224$\times$224$\times$3 |  224$\times$224$\times$3
Input size to last convolutional layer |  7$\times$7$\times$832 |  13$\times$13$\times$512 |  13$\times$13$\times$384 |  7$\times$7$\times$1024 |  14$\times$14$\times$512
Convolution filters sizes (% of conv. configs.) |  1$\times$1 (57.2%) |  1$\times$1 (71.4%) |  3$\times$3 (75%) |  1$\times$1 (66.7%) |  3$\times$3 (100%)
|  3$\times$3 (23.8%) |  3$\times$3 (28.6%) |  5$\times$5 (25%) |  3$\times$3 (33.3%) |
|  5$\times$5 (19%) |  |  |  |

Like the filter sizes, we test configurations with different input sizes and
number of filters. The input size (H and W; all of the five selected CNNs use
square inputs and filters) ranges from 7 to 224, while the number of filters
varies between 16 and 2,048. All the configurations implement a stride of 1
element and an input padding of $\frac{W_{f}-1}{2}$ and $\frac{H_{f}-1}{2}$
for the X and Y dimensions, respectively. Our experiments were performed with
data represented as single-precision floating-point values.

We compare our implementation of convolution for GPUs with those
implementations available in the NVIDIA CUDA Deep Neural Network library
(cuDNN). Among other operations used in deep neural networks, cuDNN offers
several implementations of convolution based on state–of–the–art algorithms
(GEMM, FFT, and Winograd). Most of the well-known deep learning frameworks
(like TensorFlow, PyTorch, Theano, etc.) rely on cuDNN to execute operations
on NVIDIA GPUs.

cuDNN provides several variants of each of the mentioned convolution
algorithms. See Table 2 for a detailed list with a short description of their
differences. However, since cuDNN is a closed-source library, it is difficult
to assess the exact details of the implementation of its convolutions. The
latest cuDNN version provides a helper function that uses heuristics to
suggest the fastest convolution variant for a specific convolution parameters
configuration. However, there is not much information on which heuristics are
used, and the suggested variant is not guaranteed to be the fastest. Thus, we
run all the convolution variants for all the configurations, to ensure that we
are comparing our results against the lowest cuDNN time.

Table 2: Convolution algorithm variants available in cuDNN. Algorithm | Variant | Description
---|---|---
GEMM | — | The transformed input matrix is explicitly generated before
|  | the GEMM kernel
| Implicit | The input transformation is performed on–the–fly by the kernel that
|  | computes the GEMM
| Implicit precomp. |  Like _Implicit_ , but another kernel precomputes offsets used in the
|  | implicit transformation
FFT | — | Baseline FFT-based convolution
| Tiled | The inputs are processed in tiles to reduce the temporary storage
|  | required
Winograd | — | A single kernel performs the Winograd transforms and multiplication
| Non-fused | The Winograd transform of inputs, filters and outputs is performed
|  | in separate kernels

The convolution algorithms in cuDNN experience some parameter limitations
(which are different for each algorithm), yielding some of the implementations
unavailable for certain convolution configurations. Like our implementation,
most of the cuDNN convolution variants require an additional buffer in GPU
memory to store intermediate results, with varying size depending on the
algorithm and the convolution parameters. We limit the temporary allocation
size to 1 GB. This only affects around 4% of algorithm/configuration cases,
and most of these would not be relevant even without the space limitation,
because these show poor performance compared to other algorithms for the
configurations where these require more than 1 GB of additional space.

Our test platform is an IBM POWER9 server running Red Hat Enterprise Linux
Server 7.4. The GPU installed in this system is an NVIDIA Tesla V100-SXM2
(Volta generation), and the GPU software stack available is CUDA version 9.2
with cuDNN version 7.1.

###  4.1 Speedup over cuDNN

We start by comparing the performance of our implementation with the best time
of all the convolution algorithms provided by cuDNN. For the sake of clarity,
given the number of tested configurations, we present the speedup results in
three different plots, separating the configurations by filter size (figures
5, 6 and 7). We will keep this format in the rest of this section, especially
when analyzing the reasons behind the observed performance. In fact, the
filter size is the most influential parameter and determines the best
performing cuDNN algorithm for a given configuration access . Also, even
though all configurations were tested with batch sizes up to 256, figures 5
and 6 show less batch sizes to focus on the relevant results. Note also that
the configuration parameters within a given batch size (input size, number of
filters and depth) are different in each of the plots because these depend on
the convolution configurations present in the selected CNNs. All the reported
results are the mean of nine different executions.

Figure 5: Speedup of our implementation of convolution w.r.t. the best performing cuDNN algorithm for each configuration. Configurations with $1\times 1$ filters and batch size up to 64. Labels are formatted as [inputs X&Y size]-[number of filters]-[depth]. Figure 6: Speedup of our implementation of convolution w.r.t. the best performing cuDNN algorithm for each configuration. Configurations with $3\times 3$ filters and batch size up to 16. Labels are formatted as [inputs X&Y size]-[number of filters]-[depth]. Figure 7: Speedup of our implementation of convolution w.r.t. the best performing cuDNN algorithm for each configuration. Configurations with $5\times 5$ filters and batch size up to 256. Labels are formatted as [inputs X&Y size]-[number of filters]-[depth]. Table 3: Kernel execution time ($\mu$s) for the selected configurations with 1$\times$1 filters. Algorithm | GPU kernel | A | B | C
---|---|---|---|---
GEMM |  implicit_convolve_sgemm<1024>a |  |  | 19.20
implicit |  implicit_convolve_sgemm<128>a | 128.13 | 47.87 |
| Total | 128.13 | 47.87 | 19.20
GEMM | computeOffsetsKernel | 1.98 | 2.00 | 1.89
implicit | volta_scudnn_128x64_relu_interior | 105.31 | 43.23 | 22.40
precomp. | Total | 107.29 | 45.23 | 24.29
Our | scalar_prods_kernel | 58.56 | 73.86 | 22.53
impl. | Total | 58.56 | 73.86 | 22.53
A: 7-1-1-256-832, B: 14-1-1-1024-256, C: 27-1-1-256-64.
aSome kernels show different versions depending on the configuration.

All three plots show a similar trend, where our implementation shows a clear
advantage over the best cuDNN variant for batch size of 1, although this
advantage is reduced as the batch size and the input X and Y size increase.
For configurations with 1$\times$1 filters and batch size 1, the average
speedup obtained by our implementation is 1.23$\times$, with a maximum of
2.29$\times$ for the 7-32-832 configuration (the numbers correspond to [input
size in X and Y]-[number of filters]-[depth]). Convolutions with 1$\times$1
filters are 52.3% of the tested configurations, the most common case in the
dataset.

The configurations with 3$\times$3 filters are those where our implementation
is less competitive. Except for the three leftmost configurations, the
performance attained is equal or lower than the best cuDNN variants (mainly
the two based on Winograd). The Winograd convolution was originally published
for 3$\times$3 filters Lavin15b_winograd  and obtains notably high performance
compared to the other cuDNN variants in the tested configurations with this
filter size; in around 40% of the cases the second highest performing variant
is at least 50% slower than one of the two Winograd variants.

For 5$\times$5 filters, the results are similar to the 1$\times$1 filter case.
Our implementation shows a notable advantage for batch size of 1, with an
average speedup of 1.36$\times$ and a maximum of 1.97$\times$ over the fastest
cuDNN variant. Winograd is still the fastest cuDNN algorithm for most
configurations with this filter size, although experiencing a smaller
advantage compared to the 3$\times$3 cases.

Overall, our proposed implementation is faster than the best cuDNN variant in
8.31% of the tested configurations, with an average speedup of 1.46$\times$
for these 50 configurations. Almost all of them have a batch size of 1, which
makes our implementation especially suited for the inference stage of CNN-
based applications. Given that most of the machine learning frameworks
automatically select the best-performing convolution algorithm for each
convolutional layer, our implementation will improve the performance of layers
with such configurations, without affecting the performance of the rest.

###  4.2 Algorithms Performance Comparison

To further analyze the behavior of our implementation, we profiled the
execution of a set of representative configurations to see the kernels
involved in each analyzed algorithm and their features.

Table 3 contains the kernel execution times for three selected configurations
with 1$\times$1 filters. Besides our implementation’s time, we include the two
implicit GEMMs because these are the cuDNN variants that provide the highest
performance for most of the cases with this filter size access . GEMM-impl-
precomp is the only presented variant with more than one kernel; however, the
_volta_scudnn_128x64_relu_interior_ kernel clearly dominates the total
execution time. Configuration A is one of the configurations where our
implementation shows a notable advantage over cuDNN. B (twelfth bar from the
left in Figure 5) and C are examples of configurations where the cuDNN
variants are faster than our implementation. From the profiling we observe
that our kernel is launched with a larger number of thread blocks compared to
the GEMM variants. For A, we launch 256 thread blocks, while GEMM-impl and
GEMM-impl-precomp launch 16 and 4, respectively. Thus, we expose more thread
block-level parallelism, which gives our implementation a clear advantage.
However, our approach is slower than the GEMM variants for configuration B,
where we launch 1,024 thread blocks, GEMM-impl 224 and GEMM-impl-precomp 32.
In this case, the overhead associated with thread block setup and termination
affects negatively our implementation, and hence the GEMM variants perform
higher. Configuration C features the same number of filters than A, but the
inputs size is larger, with almost 15 times more elements per input. Our
implementation’s time is close to both GEMM variants, but its performance is
more impacted by the input size increase. Since we increase the number of
threads per block to cope with a larger input size, a work-fusion optimization
at the thread level could be explored as a possible optimization (see Section
6 for further discussion).

Kernel times for two configurations with 3$\times$3 filters are presented in
Table 4. In this case, we include the two Winograd variants because those are
the highest performing cuDNN algorithm in most configurations with this filter
size. We keep the GEMM-implicit-precomp times to show the notable advantage of
Winograd even within the cuDNN convolution implementations. Compared to GEMM-
impl-precomp, Winograd non-fused is 1.93$\times$ faster for configuration A,
and Winograd is 1.67$\times$ faster for configuration B.

All the variants shown in Table 4 are composed by multiple kernels, including
our implementation, which features kernels corresponding to each of the two
stages as detailed in Section 3 (in Table 3 the second kernel does not appear
because it is not needed for configurations with 1$\times$1 filters).
Comparing the times of the two kernels, we observe that the overall time of
our implementation is dominated by the first stage. For the configurations
shown in Table 4, the second stage is 8.5% of the total time for A and only
1.14% of the total time for B.

Our implementation is clearly the fastest for configuration A. This is because
of our better exploitation of thread-block parallelism for configurations with
small inputs. Nevertheless, our implementation is also the most affected by
the increase of input size and depth from configuration B. The profiling
obtained for _scalar_prods_kernel_ shows that the compute and memory units
utilization is around 50%, even though the kernel shows an occupancy of 86%.
This further indicates an improvement direction, as we discuss later in
Section 6.

Table 4: Kernel execution time ($\mu$s) for the selected configurations with 3$\times$3 filters. Algorithm | GPU kernel | A | B
---|---|---|---
Winograd | generateWinogradTilesKernel | 9.12 | 19.77
| winograd3x3Kernel<1,4,8> | 101.91 | 212.58
| Total | 111.03 | 232.35
Winograd | winogradForwardData4x4 | 8.06 | 22.75
non-fused | winogradForwardFilter4x4 | 17.44 | 35.10
| volta_sgemm_128x64_nn | 69.31 | 242.56
| winogradForwardOutput4x4 | 10.82 | 27.14
| Total | 105.63 | 327.55
GEMM | computeOffsetsKernel | 1.98 | 2.11
implicit | volta_scudnn_128x64_relu_interior | 201.47 | 386.97
precomp. | Total | 203.45 | 389.08
Our | scalar_prods_kernel_simple | 52.86 | 461.37
impl. | sum_kernel_simple | 4.93 | 5.31
| Total | 57.79 | 466.68
A: 7-1-3-384-192, B: 13-1-3-384-384.

For configurations with 5$\times$5 filters, we selected two configurations
with different batch size to see their effect on the algorithm performance.
The execution time of the Winograd kernels in configuration A (Table 5) shows
that the pre- and post-processing kernels impacts significantly the execution
time of convolutions with relatively small computational load. In this
specific case, the sum of the _winogradForward…4x4_ kernels is 53.3% of the
total Winograd time. Even if we consider only the _volta_sgemm_128x64_nn_
time, our implementation’s time is still faster, following the trend for
1$\times$1 and 3$\times$3 filters. Comparing the times with configuration B,
we observe that Winograd scales better with the batch size than our
implementation, as it is expected from arithmetic strength reduction
optimization approaches, which are more suited for convolution configurations
with larger computational loads.

Table 5: Kernel execution time ($\mu$s) for the selected configurations with 5$\times$5 filters. Algorithm | GPU kernel | A | B
---|---|---|---
Winograd | winogradForwardData4x4 | 13.82 | 13.89
non-fused | winogradForwardFilter4x4 | 9.15 | 9.73
| volta_sgemm_128x64_nn | 34.91 | 35.36
| winogradForwardOutput4x4 | 16.92 | 17.60
| Total | 74.80 | 76.58
Our | scalar_prods_kernel | 16.80 | 107.58
impl. | sum_kernel | 5.70 | 9.02
| Total | 22.50 | 116.60
A: 7-1-5-128-48, B: 7-8-5-128-48.

##  5 Related Work

The main idea behind the GEMM-based convolution approach Chetlur_cuDNN14  is
to convert a convolution into a matrix–matrix multiplication, thus being able
to exploit already existing high-performance GEMM implementations. However,
the data transformations required in this process may be too costly for naive
implementations to be competitive (e.g. the best GEMM variant from cuDNN is
usually the implicit with offsets precomputation).

Several optimization techniques for CNNs focus on algorithmic modifications to
reduce the number of operations (mainly multiplications). As we already
introduced in Section 2, Winograd–based Lavin15b_winograd  and FFT–based
VasilacheJMCPL14_fft  convolutions are two examples of this approach. In both
cases, inputs and filter are multiplied after being transformed into another
domain, obtaining the convolution output employing less multiplication
operations. However, considering only the theoretical number of
multiplications fails to account for the possible overheads of performing the
transformations. These overheads are most noticeable for convolutions with
lower computational load, and when the transformations are implemented as
separate GPU kernels like in the non-fused variant present in cuDNN.

Another proposal within this category is to reduce the arithmetic complexity
eliminating the redundant neuron connections sparse_Han15 . For this, authors
in that paper perform an initial training to assess which weights are
irrelevant. Afterwards, they prune them by setting to zero the weights below a
certain threshold. This converts the filters into sparse volumes, which
combined with the dynamic sparsity of inputs from ReLU activations, enables
them to reduce the number of multiplications needed for the convolution.

A recent paper by Xingyu Liu et al. sparse_wino_Liu18  combines the Winograd
convolution with the pruning technique. They reorder some operations to be
able to keep inputs and filters sparse in the Winograd domain, and thus reduce
the cost of the point-wise multiplication. For the inputs, they apply the ReLU
activation function of the previous layer after the inputs are already
transformed to the Winograd domain. The filters are trained in the Winograd
domain, avoiding the need for a domain transformation. To obtain sparse
filters, they apply the pruning directly to the Winograd-domain filters.
Another work by Jongsoo Park et al. iclr_ParkLWT0CD17  also uses pruning to
improve the performance of convolution operations in CNNs. They present a
dense-matrix-with-sparse-matrix multiplication implementation for Intel Atom,
Xeon and Xeon Phy that enables the use of pruned filters in convolutional
layers, thus reducing the storage and computation cost of these layers.

A different approach is presented by Chao Li et al. mem_layout_Li16 . Instead
of focusing on the computational efficiency, they focus on the memory
efficiency of GPU implementations of convolution. As part of their work, they
study how the tensors’ data-layout affect the neural network performance, and
propose some optimizations. Similar to this work, our convolution design pays
special attention to the memory efficiency.

In a wider scope, there are several works that present other implementations
of convolution operations to improve specific scenarios of CNNs training and
deployment. For example, the work presented by Nicoli Dryden et al.
cnn_dist_train_dryden19  focuses on improving the performance of distributed
CNN training, and as part of their effort, they present a convolution
operation implementation that is better suited to exploit model-parallelism.

##  6 Conclusions and Future Work

After introducing the main features of the convolution operation used in CNNs,
we presented our design for an implementation of convolution for GPUs. As part
of this presentation, we thoroughly discussed the trade-offs involved in such
design, to finally settle with an implementation that maximizes the
exploitation of the data reuse available in the convolution operation, and
performs efficient memory accesses without the need for data-layout
transformations by selecting the most appropriate tensor layout. We evaluated
our implementation with more than 600 convolution parameter configurations
from the forward propagation phase of five of the most well known CNNs. We
compared the performance of our convolution with all of the state–of–the–art
convolution variants available in cuDNN, the reference library of GPU
operations for deep learning. The results show that our implementation is
highly competitive for certain parameter intervals, especially convolutions
with 1$\times$1 filters and small batch size, where we obtain speedups of up
to 2.29$\times$, which will contribute to the overall performance improvement
without introducing any drawbacks due to the automatic convolution algorithm
selection present in current machine learning frameworks.

As future work, we plan to explore the viability of a more sophisticated logic
to distribute the work among GPU threads and thread blocks. The execution
profiling suggests that, for some configurations, we could apply a work-fusion
optimization to improve our design. For example, fusing the work of several
thread blocks could be a potential optimization for configurations with large
number of filters and small to medium depth. However, finding the optimal work
distribution for all configurations is not straightforward, because it depends
on the specific combination of the five convolution parameters.

## References

  * (1) Abadi, M., Agarwal, A., Barham, P., Brevdo, E., Chen, Z., Citro, C., Corrado, G.S., Davis, A., Dean, J., Devin, M., Ghemawat, S., Goodfellow, I.J., Harp, A., Irving, G., Isard, M., Jia, Y., Józefowicz, R., Kaiser, L., Kudlur, M., Levenberg, J., Mané, D., Monga, R., Moore, S., Murray, D.G., Olah, C., Schuster, M., Shlens, J., Steiner, B., Sutskever, I., Talwar, K., Tucker, P.A., Vanhoucke, V., Vasudevan, V., Viégas, F.B., Vinyals, O., Warden, P., Wattenberg, M., Wicke, M., Yu, Y., Zheng, X.: TensorFlow: Large-scale machine learning on heterogeneous distributed systems.  CoRR abs/1603.04467 (2016).  URL http://arxiv.org/abs/1603.04467
  * (2) Chetlur, S., Woolley, C., Vandermersch, P., Cohen, J., Tran, J., Catanzaro, B., Shelhamer, E.: cuDNN: Efficient primitives for deep learning.  CoRR (2014)
  * (3) Dongarra, J.J., Hammarling, S., Higham, N.J., Relton, S.D., Valero-Lara, P., Zounon, M.: The design and performance of batched BLAS on modern high-performance computing systems.  In: International Conference on Computational Science (ICCS), pp. 495–504 (2017)
  * (4) Dryden, N., Maruyama, N., Moon, T., Benson, T., Snir, M., Van Essen, B.: Channel and filter parallelism for large-scale CNN training.  In: Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis, SC 2019. Association for Computing Machinery, New York, NY, USA (2019).  DOI 10.1145/3295500.3356207.  URL https://doi.org/10.1145/3295500.3356207
  * (5) Fukushima, K.: Cognitron: A self-organizing multilayer neural network.  Biological Cybernetics 20, 121–136 (1975)
  * (6) Fukushima, K., Miyake, S., Ito, T.: Neocognitron: A neural network model for a mechanism of visual pattern recognition.  IEEE Transactions on Systems, Man, and Cybernetics SMC-13, 826–834 (1983)
  * (7) Gu, J., Wang, Z., Kuen, J., Ma, L., Shahroudy, A., Shuai, B., Liu, T., Wang, X., Wang, G.: Recent advances in convolutional neural networks.  CoRR abs/1512.07108 (2015).  URL http://arxiv.org/abs/1512.07108
  * (8) Han, S., Pool, J., Tran, J., Dally, W.J.: Learning both weights and connections for efficient neural networks.  CoRR abs/1506.02626 (2015).  URL http://arxiv.org/abs/1506.02626
  * (9) He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition.  CoRR abs/1512.03385 (2015).  URL http://arxiv.org/abs/1512.03385
  * (10) Iandola, F.N., Moskewicz, M.W., Ashraf, K., Han, S., Dally, W.J., Keutzer, K.: SqueezeNet: AlexNet-level accuracy with 50x fewer parameters and <1MB model size.  CoRR abs/1602.07360 (2016).  URL http://arxiv.org/abs/1602.07360
  * (11) Jorda, M., Valero-Lara, P., Peña, A.J.: Performance evaluation of cuDNN convolution algorithms on NVIDIA Volta GPUs.  IEEE Access 7, 70461–70473 (2019).  DOI 10.1109/ACCESS.2019.2918851
  * (12) Krizhevsky, A., Sutskever, I., Hinton, G.E.: ImageNet classification with deep convolutional neural networks.  In: Advances in Neural Information Processing Systems 25, pp. 1097–1105 (2012)
  * (13) Lavin, A.: Fast algorithms for convolutional neural networks.  CoRR abs/1509.09308 (2015).  URL http://arxiv.org/abs/1509.09308
  * (14) Lecun, Y., Bottou, L., Bengio, Y., Haffner, P.: Gradient-based learning applied to document recognition.  Proceedings of the IEEE 86(11), 2278–2324 (1998)
  * (15) Li, C., Yang, Y., Feng, M., Chakradhar, S., Zhou, H.: Optimizing memory efficiency for deep convolutional neural networks on GPUs.  In: Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis (SC16), pp. 633–644 (2016).  DOI 10.1109/SC.2016.53
  * (16) Liu, X., Pool, J., Han, S., Dally, W.J.: Efficient sparse-winograd convolutional neural networks.  CoRR abs/1802.06367 (2018).  URL http://arxiv.org/abs/1802.06367
  * (17) Mathieu, M., Henaff, M., LeCun, Y.: Fast training of convolutional networks through FFTs.  CoRR abs/1312.5851 (2013).  URL http://arxiv.org/abs/1312.5851
  * (18) Park, J., Li, S.R., Wen, W., Tang, P.T.P., Li, H., Chen, Y., Dubey, P.: Faster CNNs with direct sparse convolutions and guided pruning.  In: 5th International Conference on Learning Representations, ICLR 2017, Toulon, France, April 24-26, 2017, Conference Track Proceedings (2017).  URL https://openreview.net/forum?id=rJPcZ3txx
  * (19) Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Gimelshein, N., Antiga, L., Desmaison, A., Kopf, A., Yang, E., DeVito, Z., Raison, M., Tejani, A., Chilamkurthy, S., Steiner, B., Fang, L., Bai, J., Chintala, S.: PyTorch: An imperative style, high-performance deep learning library.  In: Advances in Neural Information Processing Systems 32, pp. 8024–8035 (2019).  URL http://papers.neurips.cc/paper/9015-pytorch-an-imperative-style-high-performance-deep-learning-library.pdf
  * (20) Ramachandran, P., Zoph, B., Le, Q.V.: Searching for activation functions.  CoRR abs/1710.05941 (2017).  URL http://arxiv.org/abs/1710.05941
  * (21) Russakovsky, O., Deng, J., Su, H., Krause, J., Satheesh, S., Ma, S., Huang, Z., Karpathy, A., Khosla, A., Bernstein, M., Berg, A.C., Fei-Fei, L.: ImageNet large scale visual recognition challenge.  International Journal of Computer Vision (IJCV) 115(3), 211–252 (2015).  DOI 10.1007/s11263-015-0816-y
  * (22) Simonyan, K., Zisserman, A.: Very deep convolutional networks for large-scale image recognition.  CoRR abs/1409.1556 (2014).  URL http://arxiv.org/abs/1409.1556
  * (23) Szegedy, C., Liu, W., Jia, Y., Sermanet, P., Reed, S.E., Anguelov, D., Erhan, D., Vanhoucke, V., Rabinovich, A.: Going deeper with convolutions.  CoRR abs/1409.4842 (2014).  URL http://arxiv.org/abs/1409.4842
  * (24) Vasilache, N., Johnson, J., Mathieu, M., Chintala, S., Piantino, S., LeCun, Y.: Fast convolutional nets with fbfft: A GPU performance evaluation.  CoRR abs/1412.7580 (2014).  URL http://arxiv.org/abs/1412.7580
  * (25) Winograd, S.: Arithmetic Complexity of Computations.  CBMS-NSF Regional Conference Series in Applied Mathematics. Society for Industrial and Applied Mathematics (1980).  URL https://books.google.es/books?id=GU1NQJBcWIsC