Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-11-02 10:34:21

0001 #ifndef TMVA_SOFIE_ROPERATOR_CONV
0002 #define TMVA_SOFIE_ROPERATOR_CONV
0003 
0004 #include "TMVA/SOFIE_common.hxx"
0005 #include "TMVA/ROperator.hxx"
0006 #include "TMVA/RModel.hxx"
0007 
0008 #include <memory>
0009 #include <sstream>
0010 #include <algorithm>
0011 #include <stdexcept>
0012 #include <vector>
0013 #include <cassert>
0014 
0015 namespace TMVA {
0016 namespace Experimental {
0017 namespace SOFIE {
0018 
0019 template<typename T>
0020 class ROperator_Conv final : public ROperator
0021 {
0022 private:
0023    bool fBroadcastBias = false;
0024 
0025    std::string fAttrAutopad;
0026    std::vector<size_t> fAttrDilations;
0027    size_t fAttrGroup;
0028    std::vector<size_t> fAttrKernelShape;
0029    std::vector<size_t> fAttrPads;
0030    std::vector<size_t> fAttrStrides;
0031 
0032    std::string fNX;
0033    std::string fNW;
0034    std::string fNB;
0035    std::string fNY;
0036 
0037    std::string convK;
0038    std::string imcol;
0039 
0040    std::vector<Dim> fShapeX;
0041    std::vector<size_t> fShapeW;
0042    std::vector<size_t> fShapeB;
0043    std::vector<Dim> fShapeY;
0044 
0045    std::string fType;
0046 
0047    size_t fDim;   // dimension of the convolution
0048 
0049 
0050 public:
0051 
0052    ROperator_Conv() {}
0053 
0054    ROperator_Conv(std::string autopad, std::vector<size_t> dilations,
0055       size_t group, std::vector<size_t> kernelShape, std::vector<size_t> pads,
0056       std::vector<size_t> strides, std::string nameX, std::string nameW,
0057       std::string nameB, std::string nameY):
0058       fAttrAutopad(autopad), fAttrDilations(dilations), fAttrGroup(group), fAttrKernelShape(kernelShape),
0059       fAttrPads(pads), fAttrStrides(strides),
0060       fNX(UTILITY::Clean_name(nameX)), fNW(UTILITY::Clean_name(nameW)),
0061       fNB(UTILITY::Clean_name(nameB)), fNY(UTILITY::Clean_name(nameY))
0062    {
0063       if(std::is_same<T, float>::value) {
0064          fType = "float";
0065       } else {
0066          throw
0067             std::runtime_error("TMVA SOFIE Encountered unsupported type parsing a Conv operator");
0068       }
0069       fInputTensorNames = { fNX, fNB };
0070       fOutputTensorNames = { fNY };
0071    }
0072 
0073    ROperator_Conv(std::string autopad, std::vector<size_t> dilations,
0074       size_t group, std::vector<size_t> kernelShape, std::vector<size_t> pads,
0075       std::vector<size_t> strides, std::string nameX, std::string nameW,
0076       std::string nameY):
0077       fAttrAutopad(autopad), fAttrDilations(dilations), fAttrGroup(group), fAttrKernelShape(kernelShape),
0078       fAttrPads(pads), fAttrStrides(strides),
0079       fNX(UTILITY::Clean_name(nameX)), fNW(UTILITY::Clean_name(nameW)), fNY(UTILITY::Clean_name(nameY))
0080    {
0081       if(std::is_same<T, float>::value) {
0082          fType = "float";
0083       } else {
0084          throw
0085             std::runtime_error("TMVA SOFIE Encountered unsupported type parsing a Conv operator");
0086       }
0087       fInputTensorNames = { fNX };
0088       fOutputTensorNames = { fNY };
0089    }
0090 
0091    std::vector<ETensorType> TypeInference(std::vector<ETensorType> input) override {
0092       ETensorType out = input[0];
0093       return {out};
0094    }
0095 
0096    // function returning output shape given input
0097    std::vector<Dim> DoShapeInference(const std::vector<Dim> & input, const std::vector<size_t> & weight) {
0098       // shape of convolution input has to be (according to ONNX): N x C x H x W
0099       // Where N : batch size, C : input  channels, H : input height, W : input width
0100 
0101       if (input.size() -2 != fDim) {
0102          throw std::runtime_error("TMVA SOFIE Conv Op Shape inference - invalid input ");
0103       }
0104       if (weight.size() -2 != fDim) {
0105          throw std::runtime_error("TMVA SOFIE Conv Op Shape inference - invalid weights ");
0106       }
0107       if (fAttrGroup == 0 && input[1].isParam)
0108          throw std::runtime_error("TMVA SOFIE Conv - param shapes not supported without group attr");
0109       if (fAttrKernelShape.empty()) {
0110          if (input[2].isParam || (fDim > 1 && input[3].isParam) || (fDim > 2 && input[4].isParam))
0111             throw std::runtime_error("TMVA SOFIE Conv - param shapes not supported without kernel attr");
0112       }
0113 
0114       if (fAttrGroup == 0) {
0115          fAttrGroup = input[1].dim / weight[1];
0116       }
0117 
0118       // kernel shape
0119       size_t k1 = ((fAttrKernelShape.empty())? weight[2] : fAttrKernelShape[0]);
0120       size_t k2 = (fDim > 1) ? ((fAttrKernelShape.empty()) ? weight[3] : fAttrKernelShape[1]) : 1;
0121       size_t k3 = (fDim > 2) ? ((fAttrKernelShape.empty()) ? weight[4] : fAttrKernelShape[2]) : 1;
0122 
0123 
0124       size_t i1 = (fDim > 1) ? ((fDim > 2) ? 3 : 2) : 1;
0125       size_t i2 = (fDim > 2) ? 4 : 3;
0126       size_t i3 = 5;
0127 
0128       if (fAttrDilations.empty()) {
0129          fAttrDilations = {1, 1, 1};
0130       }
0131       fAttrDilations.resize(3);
0132       if (fDim < 3) {
0133          fAttrDilations.resize(3, 1);
0134       }
0135       // Shape of the kernel
0136       fAttrKernelShape = {k1 + (fAttrDilations[0] - 1) * (k1 - 1),
0137                           k2 + (fAttrDilations[1] - 1) * (k2 - 1),
0138                           k3 + (fAttrDilations[2] - 1) * (k3 - 1)};
0139 
0140       if (fAttrStrides.empty()) {
0141          fAttrStrides = {1, 1, 1};
0142       }
0143       if (fDim < 3)
0144          fAttrStrides.resize(3, 1);
0145 
0146       if (fAttrAutopad == "NOTSET") {
0147          if (fAttrPads.empty()) {
0148             fAttrPads = {1, 1, 1, 1, 1, 1};
0149          }
0150       } else if (fAttrAutopad == "SAME_UPPER" || fAttrAutopad == "SAME_LOWER") {
0151          for (size_t d = 0; d < fDim; ++d) {
0152             if (input[d + 2].isParam)
0153                throw std::runtime_error(
0154                   "TMVA SOFIE Conv Op: SAME padding with parametric input shape is not supported");
0155          }
0156          // ONNX SAME padding: total_pad = max(0, (ceil(in/stride)-1)*stride + kernel - in)
0157          // SAME_UPPER places extra padding at end, SAME_LOWER at beginning
0158          fAttrPads.assign(6, 0);
0159          for (size_t d = 0; d < fDim; ++d) {
0160             size_t inSize = input[d + 2].dim;
0161             size_t stride_d = fAttrStrides[d];
0162             size_t outSize = (inSize + stride_d - 1) / stride_d;
0163             int totalPad = std::max(0, (int)((outSize - 1) * stride_d + fAttrKernelShape[d]) - (int)inSize);
0164             if (fAttrAutopad == "SAME_UPPER") {
0165                fAttrPads[d] = (size_t)(totalPad / 2);
0166                fAttrPads[d + fDim] = (size_t)(totalPad - totalPad / 2);
0167             } else {
0168                fAttrPads[d] = (size_t)(totalPad - totalPad / 2);
0169                fAttrPads[d + fDim] = (size_t)(totalPad / 2);
0170             }
0171          }
0172       } else if (fAttrAutopad != "VALID") {
0173          throw
0174             std::runtime_error("TMVA SOFIE Conv Op invalid fAutopad");
0175       }
0176       // to be sure pad is vector of size 6
0177       if (fDim < 3) fAttrPads.resize(6, 0);
0178 
0179       Dim input1 = input[2];
0180       Dim input2 = (fDim > 1) ? input[3] : Dim{1};
0181       Dim input3 = (fDim > 2) ? input[4] : Dim{1};
0182 
0183       size_t pad1 = fAttrPads[0] + fAttrPads[i1];
0184 
0185       // function to get output dimension of convolution given input
0186 
0187       auto computeOutput = [&](Dim inputDim, size_t kernel, size_t pad, size_t stride) {
0188          if (!inputDim.isParam) {
0189             size_t outSize = (inputDim.dim + pad - kernel) / stride + 1;
0190             return  Dim{outSize};
0191          } else {
0192             if (stride == 1){
0193                if ((pad - kernel + 1) == 0 )
0194                   // output is same as input
0195                   return inputDim;
0196                else  {
0197                   int64_t v =  pad - kernel + 1;
0198                   std::string outStr = "(" + inputDim.param + "+" + std::to_string(v) + ")";
0199                   return Dim{ outStr, static_cast<size_t>(-1)};
0200                }
0201             } else { // general case (stride not 1)
0202                int64_t v =  pad - kernel;
0203                std::string outStr =
0204                   "((" + inputDim.param + "+" + std::to_string(v) + ")/" + std::to_string(stride) + "+1)";
0205                return Dim{ outStr, static_cast<size_t>(-1)};
0206             }
0207          }
0208          throw std::runtime_error("TMVA SOFIE Conv Op -  invalid values");
0209          return Dim{};
0210       };
0211 
0212       Dim output1 = computeOutput(input1, fAttrKernelShape[0], pad1, fAttrStrides[0]);
0213 
0214       Dim batch_size = input[0];        // first element in input tensor
0215       Dim output_channels = Dim{weight[0]};   // first element in weight tensor
0216 
0217       std::vector<Dim> ret({ batch_size, output_channels, output1 });
0218 
0219       if (fDim == 1)
0220          return ret;
0221 
0222       size_t pad2 = fAttrPads[1] + fAttrPads[i2];
0223       Dim output2 = computeOutput(input2, fAttrKernelShape[1], pad2, fAttrStrides[1]);
0224 
0225       // output is N x M x OH x OW
0226       ret.push_back(output2);
0227       if (fDim == 2)
0228          return ret;
0229 
0230       size_t pad3 = fAttrPads[2] + fAttrPads[i3];
0231       Dim output3 = computeOutput(input3, fAttrKernelShape[2], pad3, fAttrStrides[2]);
0232 
0233       // output is N x M x OH x OW x OD
0234       ret.push_back(output3);
0235       return ret;
0236    }
0237 
0238    void Initialize(RModel& model) override {
0239       fUseSession = model.UseSession();
0240       if (!model.CheckIfTensorAlreadyExist(fNX)) {
0241          throw
0242             std::runtime_error("TMVA SOFIE Conv op Input Tensor " + fNX + " is not found in model");
0243       }
0244       fShapeX = model.GetDimTensorShape(fNX);
0245       if (fShapeX.size() < 3 || fShapeX.size()  > 5) {
0246          std::cout << fNX << " : " << ConvertDimShapeToString(fShapeX) << std::endl;
0247          throw
0248             std::runtime_error("TMVA SOFIE Conv Op input data tensor" + fNX + " is not of 3,4 or 5 dimensions");
0249       }
0250       fDim = fShapeX.size() - 2;
0251       if (!model.CheckIfTensorAlreadyExist(fNW)) {
0252          throw
0253             std::runtime_error("TMVA SOFIE Conv op Input weight Tensor " + fNW + " is not found in model");
0254       }
0255       fShapeW = model.GetTensorShape(fNW);
0256       if (fShapeW.size() < 3 || fShapeW.size()  > 5) {
0257          std::cout << fNW << " : " << ConvertShapeToString(fShapeW) << std::endl;
0258          throw std::runtime_error("TMVA SOFIE Conv Op input weight tensor" + fNW + " is not of 3,4 or 5 dimensions");
0259       }
0260       fShapeY = DoShapeInference(fShapeX, fShapeW);
0261       model.AddIntermediateTensor(fNY, model.GetTensorType(fNX), fShapeY);
0262       if (fNB != "") {
0263          if (!model.CheckIfTensorAlreadyExist(fNB)) {
0264             throw
0265                std::runtime_error("TMVA SOFIE Conv op Input Tensor " + fNB + " is not found in model");
0266          }
0267          fShapeB = model.GetTensorShape(fNB);
0268          if (fShapeB.size() != 1)
0269             throw
0270                std::runtime_error("TMVA SOFIE Conv op : invalid shape for Bias tensor (is not 1D)");
0271          std::vector<Dim> targetShape(fShapeY.begin() + 1, fShapeY.end());
0272          auto shapeDimB = model.GetDimTensorShape(fNB);
0273          bool broadcast_needed = !UTILITY::AreSameShape(shapeDimB, targetShape);
0274          if (broadcast_needed) {
0275             auto original_data = model.GetInitializedTensorData(fNB);
0276             // make bias shape equal to Y shape by adding 1
0277             if (fShapeB.size() < 1)
0278                throw std::runtime_error("TMVA SOFIE Conv op: Bias Tensor has empty shape");
0279             // we assume bias tensor dimension is equal to number of filters that is the second dimension in
0280             // the output tensor
0281             if (!(shapeDimB[0] == fShapeY[1]))
0282                throw std::runtime_error("TMVA SOFIE Conv op: Bias Tensor has wrong shape: " +
0283                                            ConvertShapeToString(fShapeB));
0284             if (fType != "float")
0285                throw std::runtime_error("TMVA SOFIE Conv op: Broadcasting for non-float type tensors is not supported");
0286             // here is the actual broadcasting
0287             fBroadcastBias = true;
0288             if (!fUseSession) {
0289                // do here broadcasting
0290                std::vector<size_t> shape(fDim + 1, 1);
0291                shape[0] = fShapeB[0];
0292                auto intTargetShape = ConvertShapeToInt(targetShape);
0293                std::shared_ptr<void> new_data_ptr(
0294                   UTILITY::UnidirectionalBroadcast(static_cast<float *>(original_data.get()), shape, intTargetShape),
0295                   std::default_delete<float[]>());
0296                model.UpdateInitializedTensor(fNB, model.GetTensorType(fNB), intTargetShape, new_data_ptr);
0297                fShapeB = model.GetTensorShape(fNB);
0298             }
0299          }
0300       }
0301       // output channel size can be parametric and is an expression
0302       std::vector<Dim> outputDims = std::vector<Dim>(fShapeY.begin()+2, fShapeY.end());
0303       //check if shape is not parametric
0304       std::vector<size_t> outputInts = ConvertShapeToInt(outputDims);
0305       Dim channelDim;
0306       if (outputInts.empty()) {
0307          auto outputChannelSize = ConvertDimShapeToLength(outputDims); // size/channel = D * H * W
0308          channelDim = Dim{ outputChannelSize, static_cast<size_t>(-1)};
0309       } else {
0310          size_t outputChannelSize = ConvertShapeToLength(outputInts);
0311          channelDim = Dim{ outputChannelSize };
0312       }
0313       size_t kernelSize = fAttrKernelShape[0];
0314       for (size_t i = 1; i < fDim; i++) {
0315          kernelSize *= fAttrKernelShape[i];
0316       }
0317 
0318       std::vector<size_t> shape1 = {fShapeW[0], fShapeW[1], kernelSize};
0319       std::vector<Dim> shape2 = {Dim{fShapeW[1]}, Dim{kernelSize}, channelDim };
0320       model.AddIntermediateTensor(fNX +"_f", ConvertStringToType(fType), shape1 );
0321       model.AddIntermediateTensor(fNX +"_xcol", ConvertStringToType(fType), shape2 );
0322       convK = fNX +"_f";
0323       imcol = fNX +"_xcol";
0324       fOutputTensorNames.emplace_back(convK);
0325       fOutputTensorNames.emplace_back(imcol);
0326       fInputTensorNames.emplace_back(convK);
0327       fInputTensorNames.emplace_back(imcol);
0328 
0329       if (model.Verbose()) {
0330          std::cout << "Conv - " << fDim << "  " << fNX << " : " << ConvertDimShapeToString(fShapeX)
0331                   << " --> " << fNY << " : " << ConvertDimShapeToString(fShapeY) << std::endl;
0332       }
0333    }
0334 
0335    std::string GenerateInitCode() override {
0336       std::stringstream out;
0337       // Generate initialization code for broadcasting of bias tensor
0338       if (fBroadcastBias) {
0339          // include a separate scope to avoid defining unique operator temp variables
0340          std::vector<size_t> shape(fDim + 1, 1);
0341          // bias (is a 1D tensor)
0342          shape[0] = fShapeB[0];
0343          std::vector<Dim> targetShape(fShapeY.begin() + 1, fShapeY.end());
0344          out << "//--- broadcast bias tensor " << fNB << "for Conv op if needed \n";
0345          // in case of dynamic tensors check needs to be done at run time
0346          bool isOutDynamic = ConvertShapeToInt(targetShape).empty();
0347          auto length = ConvertDimShapeToLength(targetShape);
0348          if (isOutDynamic)
0349             out << SP << "if (" << length << " > " << ConvertShapeToLength(shape) << ") {\n";
0350          else
0351             out << SP << "{\n";
0352          out << SP << SP << "float * data = TMVA::Experimental::SOFIE::UTILITY::UnidirectionalBroadcast(tensor_"
0353              << fNB << ", " << ConvertShapeToString(shape) << ", " << ConvertDimShapeToString(fShapeY) << ");\n";
0354          out << SP << SP << "fTensor_" << fNB << ".resize(" << length << ");\n";
0355          out << SP << SP << "std::copy(data, data + " << length << ", fTensor_" << fNB << ".begin());\n";
0356          out << SP << SP << "tensor_" << fNB << " = fTensor_" << fNB << ".data();\n";
0357          out << SP << SP << "delete[] data;\n";
0358          out << SP << "}\n";
0359       }
0360       return out.str();
0361    }
0362 
0363    std::string Generate(std::string OpName) override {
0364       OpName = "op_" + OpName;
0365 
0366       if (fShapeX.empty() || fShapeW.empty() || (fNB != "" && fShapeB.empty()) || fShapeY.empty()) {
0367          throw
0368             std::runtime_error("TMVA SOFIE Conv Op called to Generate without being initialized first");
0369       }
0370 
0371       std::stringstream out;
0372       auto bsize = fShapeX[0];
0373       size_t kDepth = (fDim > 2) ?  fShapeW[2] : 1;  // kernel depth
0374       size_t kHeight = (fDim > 1) ? fShapeW[fDim] : 1;  // kernel height
0375       size_t kWidth = fShapeW[fDim+1]; // kernel width
0376       auto iDepth = (fDim > 2) ?  fShapeX[2] : Dim{1};  // input depth
0377       auto iHeight = (fDim > 1) ? fShapeX[fDim] : Dim{1}; // input height
0378       auto iWidth = fShapeX[fDim+1]; // input width
0379       auto oDepth = (fDim > 2) ? fShapeY[2] : Dim{1}; // output depth
0380       auto oHeight = (fDim > 1) ? fShapeY[fDim] : Dim{1};  // ouput height
0381       auto oWidth = fShapeY[fDim+1]; // output width
0382       // total output size for a channel
0383       auto outputChannelStride = ConvertDimShapeToLength(std::vector<Dim>{oDepth, oHeight, oWidth}); // size of channel = D * H * W
0384       auto outputBatchStride =  ConvertDimShapeToLength(std::vector<Dim>{fShapeY[1] , oDepth, oHeight, oWidth}); // size of C * D * H * W
0385       // input size
0386       auto inputChannelStride = ConvertDimShapeToLength(std::vector<Dim>{iDepth, iHeight, iWidth});
0387       auto inputBatchStride =  ConvertDimShapeToLength(std::vector<Dim>{fShapeX[1] , iDepth, iHeight, iWidth}); // size of C * D * H * W
0388 
0389       out << "\n//----  operator Conv " << OpName << "\n";
0390 
0391       // vectorize the (dilated)convolution kernels into a matrix
0392       // no need to transpose the matrix
0393       // to fix for 1d and 3d
0394 
0395       size_t id = (fDim > 2) ? fDim-3 : 2;
0396       size_t ih = (fDim > 1) ? fDim-2 : 1;
0397       size_t iw = fDim-1;
0398 
0399       size_t wstrideDil = fAttrDilations[iw];
0400       size_t hstride = kWidth;
0401       size_t hstrideDil = fAttrDilations[ih] * fAttrKernelShape[iw];  // stride dilated in the height
0402       size_t dstride = kHeight * kWidth;
0403       size_t dstrideDil = fAttrDilations[id] * fAttrKernelShape[ih] * fAttrKernelShape[iw];
0404       size_t icstride = kHeight * kWidth * kDepth;
0405       size_t icstrideDil = fAttrKernelShape[id] * fAttrKernelShape[ih] * fAttrKernelShape[iw];
0406       size_t ocstride = fShapeW[1] * icstride;
0407       size_t ocstrideDil = fShapeW[1] * icstrideDil;
0408 
0409       out << SP << "for (std::size_t oc = 0; oc < " << fShapeW[0] << "; oc++) {\n";
0410       out << SP << SP << "for (std::size_t ic = 0; ic < " << fShapeW[1] << "; ic++) {\n";
0411       if (fDim > 2)
0412          out << SP << SP << SP << "for (std::size_t kd = 0; kd < " << kDepth << "; kd++) {\n";
0413       if (fDim > 1)
0414          out << SP << SP << SP << "for (std::size_t kh = 0; kh < " << kHeight << "; kh++) {\n";
0415       out << SP << SP << SP << SP << "for (std::size_t kw = 0; kw < " << kWidth << "; kw++) {\n";
0416 
0417       out << SP << SP << SP << SP << SP << "tensor_" <<fNX <<  "_f[oc * "
0418           << ocstrideDil << " + ic * " << icstrideDil;
0419       if (fDim > 2) out << " + kd * " << dstrideDil;
0420       if (fDim > 1) out << " + kh * " << hstrideDil;
0421       out << " + kw * " << wstrideDil  << "  ] = tensor_" << fNW << "[oc * " << ocstride << " + ic * " << icstride;
0422       if (fDim > 2) out << " + kd * " << dstride;
0423       if (fDim > 1) out << " + kh * " << hstride;
0424       out  << " + kw ];\n";
0425 
0426       out << SP << SP << SP << SP << "}\n";
0427       if (fDim > 1) out << SP << SP << SP << "}\n";
0428       if (fDim > 2) out << SP << SP << SP << "}\n";
0429       out << SP << SP << "}\n";
0430       out << SP << "}\n";
0431 
0432       // Dilation is already folded into the expanded kernel shape and the dilated tensor_<X>_f
0433       // layout above, so the dense im2col below must use dilation 1 to avoid double-counting it.
0434       fAttrDilations = std::vector<size_t>(3, 1);
0435 
0436       //out << SP << "char " << OpName << "_transA = 'T';\n";
0437       out << SP << "char " << OpName << "_transA = 'N';\n";
0438       out << SP << "char " << OpName << "_transB = 'N';\n";
0439       out << SP << "int " << OpName << "_m = " << outputChannelStride << ";\n"; // output h*w
0440       assert(fShapeY[1] == fShapeW[0]);
0441       //assert(fShapeW[1] == fShapeX[1] / fAttrGroup);
0442       out << SP << "int " << OpName << "_n = " << fShapeW[0] << ";\n"; // output channels
0443       out << SP << "int " << OpName << "_k = " << fShapeW[1] * fAttrKernelShape[0] * fAttrKernelShape[1] * fAttrKernelShape[2] << ";\n";
0444       out << SP << "float " << OpName << "_alpha = 1.0;\n";
0445       if (fNB != "")
0446          out << SP << "float " << OpName << "_beta = 1.0;\n";
0447       else // when bias is not present beta needs to be equal to zero to avoid re-using previous results in output tensor
0448          out << SP << "float " << OpName << "_beta = 0.0;\n";
0449 
0450 
0451       // Loop on batch size
0452       out << SP << "for (size_t n = 0; n < " << bsize << "; n++) {\n";
0453 
0454       // IM2COL: Unroll the input tensor
0455       // order input data as  (e.g. kernel 2x2)  and (xa,ya) is channel 1 and (xb,yb) is channel 2
0456       //   (xa1,..,xak,ya1,..yak)(xb1,...,xbk,yb1,..,ybk)
0457       //   (xa2,...xak+1,ya1,...yak)(......)
0458       // trick for speed is using caffe im2col and output a matrix which contains filtered values as rows.
0459       // By doing this one has consecutive memory reads and writes
0460       // Resulting matrix op_xcol is (input channels * filter_h * filter_w , output_h * output_w)
0461       if (fDim ==1) {
0462          if (fAttrPads[0] != fAttrPads[1] ) {
0463             std::cout << "TMVA SOFIE Operator Conv:  asymmetric padding not supported. Assume an average padding "
0464                       << std::endl;
0465             fAttrPads[0] = (fAttrPads[0] + fAttrPads[1]) / 2;
0466          }
0467          fAttrPads[1] = 0;
0468          fAttrStrides[1] = 1;
0469       }
0470       if (fDim == 2) {
0471          if (fAttrPads[0] != fAttrPads[2] || fAttrPads[1] != fAttrPads[3]) {
0472             std::cout << "TMVA SOFIE Operator Conv:  asymmetric padding not supported. Assume an average padding " << std::endl;
0473             fAttrPads[0] = (fAttrPads[0] + fAttrPads[2]) / 2;
0474             fAttrPads[1] = (fAttrPads[1] + fAttrPads[3]) / 2;
0475          }
0476       }
0477       if (fDim == 3) {
0478          if (fAttrPads[0] != fAttrPads[3] || fAttrPads[1] != fAttrPads[4] || fAttrPads[2] != fAttrPads[5]) {
0479             std::cout << "TMVA SOFIE Operator Conv:  asymmetric padding not supported. Assume an average padding " << std::endl;
0480             fAttrPads[0] = (fAttrPads[0] + fAttrPads[3]) / 2;
0481             fAttrPads[1] = (fAttrPads[1] + fAttrPads[4]) / 2;
0482             fAttrPads[2] = (fAttrPads[2] + fAttrPads[5]) / 2;
0483          }
0484       }
0485       out << SP << SP << "size_t out_offset = n * " << outputBatchStride  << ";\n";
0486 
0487       if (fAttrGroup == 1) {
0488          out << SP << SP << "size_t x_offset = n * " << inputBatchStride << ";\n";
0489          // when using im2col - resulting matrix is transposed, the dimension is (input_c * filter_h * filter_y,  output_h *
0490          // output_w)
0491          if (fDim < 3) {
0492             out << SP << SP << "TMVA::Experimental::SOFIE::UTILITY::Im2col<float>(tensor_" << fNX
0493                 << " + x_offset,"
0494                 //  channels, height, width, kernel_h, kernel_w, pad_h, pad_w, stride_h, stride_w, dilation_h,
0495                 //  dilation_w,
0496                 //
0497                 << fShapeW[1] << "," << iHeight << "," << iWidth << ",";
0498             if (fDim == 1)
0499                out << "1, " << fAttrKernelShape[0] << ",0," << fAttrPads[0] << ",1," << fAttrStrides[0] << ",1,"
0500                    << fAttrDilations[0];
0501             else // dim ==2
0502                out << fAttrKernelShape[0] << "," << fAttrKernelShape[1] << "," << fAttrPads[0] << "," << fAttrPads[1]
0503                    << "," << fAttrStrides[0] << "," << fAttrStrides[1] << "," << fAttrDilations[0] << ","
0504                    << fAttrDilations[1];
0505             out << "," << "tensor_" <<fNX << "_xcol);\n\n ";
0506          } else {
0507             // 3d im2col
0508             out << SP << SP << "TMVA::Experimental::SOFIE::UTILITY::Im2col_3d<float>(tensor_" << fNX
0509                 << " + x_offset,"
0510                 //  channels, d, h, w, k_d, k_h, k_w, pad_d, pad_h, pad_w, stride_d, stride_h, stride_w,
0511                 //  dilation_d, dilation_h, dilation_w,
0512                 //
0513                 << fShapeW[1] << "," << iDepth << "," << iHeight << "," << iWidth << ","
0514                 << fAttrKernelShape[0] << "," << fAttrKernelShape[1] << "," << fAttrKernelShape[2] << ","
0515                 << fAttrPads[0] << "," << fAttrPads[1] << "," << fAttrPads[2] << ","
0516                 << fAttrStrides[0] << "," << fAttrStrides[1] << "," << fAttrStrides[2] << ","
0517                 << fAttrDilations[0] << "," << fAttrDilations[1] << "," << fAttrDilations[2] << ","
0518                 << "tensor_" << fNX << "_xcol);\n\n ";
0519          }
0520          // BLAS
0521          out << SP << "TMVA::Experimental::SOFIE::Gemm_Call("
0522              << "tensor_" << fNY << " + out_offset, false, false, "
0523              << OpName << "_m, " << OpName << "_n, " << OpName << "_k, "
0524              << OpName << "_alpha, " << "tensor_" << fNX << "_xcol, tensor_" << fNX << "_f, "
0525              << OpName << "_beta, ";
0526          if (fNB != "")
0527             out << "tensor_" << fNB;
0528          else
0529             out << "nullptr";
0530          out << ");\n";
0531 
0532 
0533          // out << SP << SP << "BLAS::sgemm_(&" << OpName << "_transA, &" << OpName << "_transB, &" << OpName << "_m, &"
0534          //     << OpName << "_n, &" << OpName << "_k, &" << OpName << "_alpha, " << "tensor_" << fNX << "_xcol, &" << OpName
0535          //     << "_m,\n"; // use m if op_xcol is not transpose , otherwise k
0536          // out << SP << SP << SP << "tensor_" << fNX << "_f, &" << OpName << "_k, &" << OpName << "_beta, tensor_" << fNY
0537          //     << " + out_offset, &" << OpName << "_m);\n";
0538       } else {
0539          // case of group convolution
0540          // Unroll (IM2COL) the input tensor- make loop on groups and repeat operations (IM2COL + GEMM for each
0541          // group)
0542          // out << SP << SP << "size_t out_offset = n * " << fShapeY[1] * oDepth * oHeight * oWidth << ";\n";
0543          out << SP << SP << "for (size_t g = 0; g < " << fAttrGroup << "; g++) {\n";
0544          out << SP << SP << "size_t x_offset = n * " << inputBatchStride << " + g * "
0545              << fShapeW[1] << " * " << inputChannelStride << ";\n ";
0546          out << SP << SP << "size_t g_offset = g * " << fShapeW[0] << " * (" << outputChannelStride << ") / " << fAttrGroup << ";\n ";
0547          out << SP << SP << "size_t out_offset = n * " << outputBatchStride << " + g_offset;\n";
0548 
0549          if (fDim < 3) {
0550             out << SP << SP << "TMVA::Experimental::SOFIE::UTILITY::Im2col<float>(tensor_" << fNX
0551                 << " + x_offset,"
0552                 //  channels, height, width, kernel_h, kernel_w, pad_h, pad_w, stride_h, stride_w, dilation_h,
0553                 //  dilation_w,
0554                 //
0555                 << fShapeW[1] << "," << iHeight << "," << iWidth << ",";
0556             if (fDim == 1)
0557                out << "1, " << fAttrKernelShape[0] << ",0," << fAttrPads[0] << ",1," << fAttrStrides[0] << ",1,"
0558                    << fAttrDilations[0];
0559             else // dim ==2
0560                out << fAttrKernelShape[0] << "," << fAttrKernelShape[1] << "," << fAttrPads[0] << "," << fAttrPads[1]
0561                    << "," << fAttrStrides[0] << "," << fAttrStrides[1] << "," << fAttrDilations[0] << ","
0562                    << fAttrDilations[1];
0563             out << ", tensor_" << fNX << "_xcol);\n\n ";
0564          } else {
0565             // 3d im2col
0566             out << SP << SP << "TMVA::Experimental::SOFIE::UTILITY::Im2col_3d<float>(tensor_" << fNX
0567                 << " + x_offset,"
0568                 //  channels, d, h, w, k_d, k_h, k_w, pad_d, pad_h, pad_w, stride_d, stride_h, stride_w,
0569                 //  dilation_d, dilation_h, dilation_w,
0570                 //
0571                 << fShapeW[1] << "," << iDepth << "," << iHeight << "," << iWidth << "," << fAttrKernelShape[0] << ","
0572                 << fAttrKernelShape[1] << "," << fAttrKernelShape[2] << "," << fAttrPads[0] << "," << fAttrPads[1]
0573                 << "," << fAttrPads[2] << "," << fAttrStrides[0] << "," << fAttrStrides[1] << "," << fAttrStrides[2]
0574                 << "," << fAttrDilations[0] << "," << fAttrDilations[1] << "," << fAttrDilations[2] << ",tensor_" << fNX
0575                 << "_xcol);\n\n ";
0576          }
0577 
0578          // BLAS
0579          // n must be divided by the number of groups
0580          out << SP << SP << SP << OpName << "_n = " << fShapeW[0] / fAttrGroup << ";\n";
0581          // offset g must be  g * k * n
0582          out << SP << SP << SP << "size_t offset_f = g * "
0583              << fShapeW[0] * fShapeW[1] * fAttrKernelShape[0] * fAttrKernelShape[1] * fAttrKernelShape[2] / fAttrGroup
0584              << ";\n";
0585 
0586          out << SP << "TMVA::Experimental::SOFIE::Gemm_Call("
0587              << "tensor_" << fNY << " + out_offset, false, false, "
0588              << OpName << "_m, " << OpName << "_n, " << OpName << "_k, "
0589              << OpName << "_alpha, " << "tensor_" << fNX << "_xcol, tensor_" << fNX << "_f + offset_f, "
0590              << OpName << "_beta, ";
0591          if (fNB != "")
0592             out << "tensor_" << fNB << " + g_offset";
0593          else
0594             out << "nullptr";
0595          out << ");\n";
0596 
0597          // out << SP << SP << "BLAS::sgemm_(&" << OpName << "_transA, &" << OpName << "_transB, &" << OpName << "_m, &"
0598          //     << OpName << "_n, &" << OpName << "_k, &" << OpName << "_alpha, tensor_" << fNX << "_xcol, &" << OpName
0599          //     << "_m,\n"; // use m if op_xcol is not transpose , otherwise k
0600          // out << SP << SP << SP << "tensor_" << fNX << "_f + offset_f, &" << OpName << "_k, &" << OpName << "_beta, tensor_" << fNY
0601          //     << " + out_offset"
0602          //     << ", &" << OpName << "_m);\n";
0603 
0604          out << SP << SP << "}\n"; // end of group loop
0605       }
0606 
0607       // if (fNB != "") {
0608       //    out << SP << "int " << OpName << "_size = " << outputBatchStride << ";\n";
0609       //    out << SP << "float " << OpName << "_gamma = 1.0;\n";
0610       //    out << SP << "int " << OpName << "_incx = 1;\n";
0611       //    out << SP << "int " << OpName << "_incy = 1;\n";
0612 
0613       //    out << SP << "BLAS::saxpy_(&" << OpName << "_size, &" << OpName << "_gamma, tensor_" << fNB << ", &"
0614       //        << OpName << "_incx, tensor_" << fNY << " + out_offset, &" << OpName << "_incy);\n";
0615 
0616       // }
0617       out << SP << "}\n"; // end of batch size loop
0618 
0619       return out.str();
0620       }
0621 
0622    /*! \brief Returns the blas routines needed to compile the generated code
0623     */
0624    std::vector<std::string> GetBlasRoutines() override { return { std::string("Gemm"), std::string("Axpy") }; }
0625 };
0626 
0627 } // namespace SOFIE
0628 } // namespace Experimental
0629 } // namespace TMVA
0630 
0631 #endif