1 /*
2  * Copyright (c) 2022-2023 Arm Limited.
3  *
4  * SPDX-License-Identifier: MIT
5  *
6  * Permission is hereby granted, free of charge, to any person obtaining a copy
7  * of this software and associated documentation files (the "Software"), to
8  * deal in the Software without restriction, including without limitation the
9  * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
10  * sell copies of the Software, and to permit persons to whom the Software is
11  * furnished to do so, subject to the following conditions:
12  *
13  * The above copyright notice and this permission notice shall be included in all
14  * copies or substantial portions of the Software.
15  *
16  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
17  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
18  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
19  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
20  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
21  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
22  * SOFTWARE.
23  */
24 
25 #include <algorithm>
26 #include <cstddef>
27 
28 #include <arm_neon.h>
29 
30 namespace arm_conv {
31 namespace winograd {
32 namespace output_transform {
33 
arm_fp32_1x6_1x3(unsigned int n_channels,const float * inptr,size_t matrix_stride,const float * bptr,float * outptr,size_t,size_t output_col_stride,float output_min,float output_max)34 void arm_fp32_1x6_1x3(
35   unsigned int n_channels,
36   const float* inptr,
37   size_t matrix_stride,
38   const float* bptr,
39   float *outptr,
40   size_t,  // No need to stride across rows
41   size_t output_col_stride,
42   float output_min,
43   float output_max
44 )
45 {
46   constexpr unsigned int inner_tile_cols = 8, output_tile_cols = 6;
47 
48   // For each channel of the output
49   for (; n_channels >= 4; n_channels -= 4)
50   {
51     // Matrices used and computed during this transform
52     float32x4_t F[inner_tile_cols], f[output_tile_cols], b = vdupq_n_f32(0.0f);
53 
54     // Read a 1x8 tile in the Winograd domain
55     for (auto j = 0u; j < inner_tile_cols; j++)
56     {
57       F[j] = vld1q_f32(inptr + j*matrix_stride);
58     }
59     inptr += 4;
60 
61     f[0] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[6], 1), F[5], 1), F[4], 1), F[3], 1), F[2], 1), F[1], 1), F[0], 1);
62     f[1] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[2], 1), F[6], 3), F[4], 2), F[3], -2), F[5], -3), F[1], -1);
63     f[2] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[2], 1), F[1], 1), F[6], 9), F[5], 9), F[4], 4), F[3], 4);
64     f[3] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[2], 1), F[6], 27), F[4], 8), F[3], -8), F[5], -27), F[1], -1);
65     f[4] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[2], 1), F[1], 1), F[6], 81), F[5], 81), F[4], 16), F[3], 16);
66     f[5] = vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmlaq_n_f32(vmulq_n_f32(F[7], 1), F[2], 1), F[6], 243), F[4], 32), F[3], -32), F[5], -243), F[1], -1);
67 
68     // Write out the output tile
69     if (bptr != 0)
70     {
71       b = vld1q_f32(bptr);
72       bptr += 4;
73     }
74     for (auto j = 0u; j < output_tile_cols; j++)
75     {
76       const auto y = vminq_f32(vmaxq_f32(f[j] + b, vdupq_n_f32(output_min)),
77                                vdupq_n_f32(output_max));
78       vst1q_f32(outptr + j*output_col_stride, y);
79     }
80     outptr += 4;
81   }
82   for (; n_channels >= 2; n_channels -= 2)
83   {
84     // Matrices used and computed during this transform
85     float32x2_t F[inner_tile_cols], f[output_tile_cols], b = vdup_n_f32(0.0f);
86 
87     // Read a 1x8 tile in the Winograd domain
88     for (auto j = 0u; j < inner_tile_cols; j++)
89     {
90       F[j] = vld1_f32(inptr + j*matrix_stride);
91     }
92     inptr += 2;
93 
94     f[0] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[6], 1), F[5], 1), F[4], 1), F[3], 1), F[2], 1), F[1], 1), F[0], 1);
95     f[1] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[2], 1), F[6], 3), F[4], 2), F[3], -2), F[5], -3), F[1], -1);
96     f[2] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[2], 1), F[1], 1), F[6], 9), F[5], 9), F[4], 4), F[3], 4);
97     f[3] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[2], 1), F[6], 27), F[4], 8), F[3], -8), F[5], -27), F[1], -1);
98     f[4] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[2], 1), F[1], 1), F[6], 81), F[5], 81), F[4], 16), F[3], 16);
99     f[5] = vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmla_n_f32(vmul_n_f32(F[7], 1), F[2], 1), F[6], 243), F[4], 32), F[3], -32), F[5], -243), F[1], -1);
100 
101     // Write out the output tile
102     if (bptr != 0)
103     {
104       b = vld1_f32(bptr);
105       bptr += 2;
106     }
107     for (auto j = 0u; j < output_tile_cols; j++)
108     {
109       const auto y = vmin_f32(vmax_f32(f[j] + b, vdup_n_f32(output_min)),
110                               vdup_n_f32(output_max));
111       vst1_f32(outptr + j*output_col_stride, y);
112     }
113     outptr += 2;
114   }
115   for (; n_channels; n_channels--)
116   {
117     // Matrices used and computed during this transform
118     float F[inner_tile_cols], f[output_tile_cols], b = 0.0f;
119 
120     // Read a 1x8 tile in the Winograd domain
121     for (auto j = 0u; j < inner_tile_cols; j++)
122     {
123       F[j] = *(inptr + j*matrix_stride);
124     }
125     inptr++;
126 
127     f[0] = F[0]*1 + F[1]*1 + F[2]*1 + F[3]*1 + F[4]*1 + F[5]*1 + F[6]*1;
128     f[1] = F[1]*-1 + F[5]*-3 + F[3]*-2 + F[4]*2 + F[6]*3 + F[2]*1;
129     f[2] = F[3]*4 + F[4]*4 + F[5]*9 + F[6]*9 + F[1]*1 + F[2]*1;
130     f[3] = F[1]*-1 + F[5]*-27 + F[3]*-8 + F[4]*8 + F[6]*27 + F[2]*1;
131     f[4] = F[3]*16 + F[4]*16 + F[5]*81 + F[6]*81 + F[1]*1 + F[2]*1;
132     f[5] = F[1]*-1 + F[5]*-243 + F[3]*-32 + F[4]*32 + F[6]*243 + F[2]*1 + F[7]*1;
133 
134     // Write out the output tile
135     if (bptr != 0)
136     {
137       b = *(bptr++);
138     }
139     for (auto j = 0u; j < output_tile_cols; j++)
140     {
141       *(outptr + j*output_col_stride) = std::max(std::min(f[j] + b, output_max), output_min);
142     }
143     outptr++;
144   }
145 }
146 
147 }  // namespace output_transform
148 }  // namespace winograd
149 }  // namespace arm_conv
150