1 //------------------------------------------------------------------------------
2 // LAGraph_dnn: sparse deep neural network
3 //------------------------------------------------------------------------------
4
5 /*
6 LAGraph: graph algorithms based on GraphBLAS
7
8 Copyright 2020 LAGraph Contributors.
9
10 (see Contributors.txt for a full list of Contributors; see
11 ContributionInstructions.txt for information on how you can Contribute to
12 this project).
13
14 All Rights Reserved.
15
16 NO WARRANTY. THIS MATERIAL IS FURNISHED ON AN "AS-IS" BASIS. THE LAGRAPH
17 CONTRIBUTORS MAKE NO WARRANTIES OF ANY KIND, EITHER EXPRESSED OR IMPLIED,
18 AS TO ANY MATTER INCLUDING, BUT NOT LIMITED TO, WARRANTY OF FITNESS FOR
19 PURPOSE OR MERCHANTABILITY, EXCLUSIVITY, OR RESULTS OBTAINED FROM USE OF
20 THE MATERIAL. THE CONTRIBUTORS DO NOT MAKE ANY WARRANTY OF ANY KIND WITH
21 RESPECT TO FREEDOM FROM PATENT, TRADEMARK, OR COPYRIGHT INFRINGEMENT.
22
23 Released under a BSD license, please see the LICENSE file distributed with
24 this Software or contact permission@sei.cmu.edu for full terms.
25
26 Created, in part, with funding and support from the United States
27 Government. (see Acknowledgments.txt file).
28
29 This program includes and/or can make use of certain third party source
30 code, object code, documentation and other files ("Third Party Software").
31 See LICENSE file for more details.
32
33 */
34
35 //------------------------------------------------------------------------------
36
37 // LAGraph_dnn: sparse deep neural network. Contributed by Tim Davis,
38 // Texas A&M University. Based on inferenceReLUvec.m by Jeremy Kepner, MIT.
39
40 // Performs ReLU inference using input feature vectors Y0.
41
42 // See http://graphchallenge.org/ for a description of the algorithm.
43
44 // On input, Y0 is the initial feature vectors, of size nfeatures-by-nneurons.
45 // This format uses the graph convention that A(i,j) is the edge (i,j),
46 // Each row of Y0 is a single feature.
47
48 // W is an array of size nlayers of sparse matrices. Each W[layer] matrix has
49 // the same size: nneurons-by-nneurons. W[layer] represents the DNN weights
50 // for that layer.
51
52 // The Bias[layer] matrices are diagaonal, and the same size as W[layer].
53
54 // All matrices must have the same type: either GrB_FP32 or GrB_FP64.
55
56 // On output, Y is the computed result, of the same size and type as Y0.
57
58 #include "LAGraph.h"
59
60 #define LAGRAPH_FREE_ALL \
61 { \
62 GrB_free (&M) ; \
63 GrB_free (Yhandle) ; \
64 }
65
LAGraph_dnn(GrB_Matrix * Yhandle,GrB_Matrix * W,GrB_Matrix * Bias,int nlayers,GrB_Matrix Y0)66 GrB_Info LAGraph_dnn // returns GrB_SUCCESS if successful
67 (
68 // output
69 GrB_Matrix *Yhandle, // Y, created on output
70 // input: not modified
71 GrB_Matrix *W, // W [0..nlayers-1], each nneurons-by-nneurons
72 GrB_Matrix *Bias, // Bias [0..nlayers-1], diagonal nneurons-by-nneurons
73 int nlayers, // # of layers
74 GrB_Matrix Y0 // input features: nfeatures-by-nneurons
75 )
76 {
77
78 //--------------------------------------------------------------------------
79 // check inputs
80 //--------------------------------------------------------------------------
81
82 GrB_Info info ;
83 if (Yhandle == NULL || W == NULL || Bias == NULL || Y0 == NULL)
84 {
85 return (GrB_NULL_POINTER) ;
86 }
87
88 GrB_Matrix Y = NULL ;
89 GrB_Matrix M = NULL ;
90 (*Yhandle) = NULL ;
91
92 GrB_Type type ;
93 LAGRAPH_OK (GxB_Matrix_type (&type, Y0)) ;
94
95 GrB_Semiring plus_times, plus_plus ;
96 GrB_UnaryOp gt0, id, ymax ;
97
98 if (type == GrB_FP32)
99 {
100 plus_times = LAGraph_PLUS_TIMES_FP32 ;
101 plus_plus = LAGraph_PLUS_PLUS_FP32 ;
102 gt0 = LAGraph_GT0_FP32 ;
103 ymax = LAGraph_YMAX_FP32 ;
104 id = GrB_IDENTITY_FP32 ;
105 }
106 else if (type == GrB_FP64)
107 {
108 plus_times = LAGraph_PLUS_TIMES_FP64 ;
109 plus_plus = LAGraph_PLUS_PLUS_FP64 ;
110 gt0 = LAGraph_GT0_FP64 ;
111 ymax = LAGraph_YMAX_FP64 ;
112 id = GrB_IDENTITY_FP64 ;
113 }
114 else
115 {
116 return (GrB_DOMAIN_MISMATCH) ;
117 }
118
119 for (int layer = 0 ; layer < nlayers ; layer++)
120 {
121 GrB_Type type2 ;
122 LAGRAPH_OK (GxB_Matrix_type (&type2, W [layer])) ;
123 if (type != type2) return (GrB_DOMAIN_MISMATCH) ;
124 LAGRAPH_OK (GxB_Matrix_type (&type2, Bias [layer])) ;
125 if (type != type2) return (GrB_DOMAIN_MISMATCH) ;
126 }
127
128 //--------------------------------------------------------------------------
129 // create the output matrix Y
130 //--------------------------------------------------------------------------
131
132 GrB_Index nfeatures, nneurons ;
133 LAGRAPH_OK (GrB_Matrix_nrows (&nfeatures, Y0)) ;
134 LAGRAPH_OK (GrB_Matrix_ncols (&nneurons, Y0)) ;
135 LAGRAPH_OK (GrB_Matrix_new (&Y, type, nfeatures, nneurons)) ;
136 LAGRAPH_OK (GrB_Matrix_new (&M, GrB_BOOL, nfeatures, nneurons)) ;
137
138 //--------------------------------------------------------------------------
139 // propagate the features through the neuron layers
140 //--------------------------------------------------------------------------
141
142 // double t1 = 0, t2 = 0, t3 = 0, t4 = 0, t ;
143
144 for (int layer = 0 ; layer < nlayers ; layer++)
145 {
146 // Y = Y * W [layer], using the conventional PLUS_TIMES semiring
147 // t = omp_get_wtime ( ) ;
148 LAGRAPH_OK (GrB_mxm (Y, NULL, NULL, plus_times,
149 ((layer == 0) ? Y0 : Y), W [layer], NULL)) ;
150 // t1 += (omp_get_wtime ( ) - t) ;
151
152 // Y = Y * Bias [layer], using the PLUS_PLUS semiring. This computes
153 // Y(i,j) += Bias [layer] (j,j) for each entry Y(i,j). It does not
154 // introduce any new entries in Y.
155 // t = omp_get_wtime ( ) ;
156 LAGRAPH_OK (GrB_mxm (Y, NULL, NULL, plus_plus, Y, Bias [layer], NULL)) ;
157 // t2 += (omp_get_wtime ( ) - t) ;
158
159 // delete entries from Y: keep only those entries greater than zero
160 #if defined ( GxB_SUITESPARSE_GRAPHBLAS ) \
161 && ( GxB_IMPLEMENTATION >= GxB_VERSION (3,0,0) )
162 // using SuiteSparse:GraphBLAS 3.0.0 or later.
163 // t = omp_get_wtime ( ) ;
164 LAGRAPH_OK (GxB_select (Y, NULL, NULL, GxB_GT_ZERO, Y, NULL, NULL)) ;
165 // t3 += (omp_get_wtime ( ) - t) ;
166
167 #else
168 // using SuiteSparse v2.x or earlier, or any other GraphBLAS library.
169 LAGRAPH_OK (GrB_apply (M, NULL, NULL, gt0, Y, NULL)) ;
170 LAGRAPH_OK (GrB_apply (Y, M, NULL, id, Y, LAGraph_desc_ooor)) ;
171 #endif
172
173 // threshold maximum values: Y (Y > 32) = 32
174 // t = omp_get_wtime ( ) ;
175 LAGRAPH_OK (GrB_apply (Y, NULL, NULL, ymax, Y, NULL)) ;
176 // t4 += (omp_get_wtime ( ) - t) ;
177 }
178
179 // printf ("\nY*W: %g sec\n", t1) ;
180 // printf ("Y+B: %g sec\n", t2) ;
181 // printf ("RelU %g sec\n", t3) ;
182 // printf ("Ymax %g sec\n", t4) ;
183
184 //--------------------------------------------------------------------------
185 // free workspace and return result
186 //--------------------------------------------------------------------------
187
188 GrB_free (&M) ;
189 (*Yhandle) = Y ;
190 return (GrB_SUCCESS) ;
191 }
192
193