[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r314224603
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,302 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (DType), size of which is determined by a workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report
+// `Auto-Differentiating Linear Algebra` for details
+// on https://arxiv.org/pdf/1710.08717.pdf
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+// dA overwritten by L^-1 dA
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+// X (square) overwritten by X L
+// Y overwritten by the diagonal of X
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313780388
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313822055
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,302 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (DType), size of which is determined by a workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report
+// `Auto-Differentiating Linear Algebra` for details
+// on https://arxiv.org/pdf/1710.08717.pdf
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+// dA overwritten by L^-1 dA
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+// X (square) overwritten by X L
+// Y overwritten by the diagonal of X
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313822055
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,302 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (DType), size of which is determined by a workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report
+// `Auto-Differentiating Linear Algebra` for details
+// on https://arxiv.org/pdf/1710.08717.pdf
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+// dA overwritten by L^-1 dA
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+// X (square) overwritten by X L
+// Y overwritten by the diagonal of X
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313816473
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
Review comment:
Yes. Cited.
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313782597
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313782129
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
Review comment:
I tried it out. But I found std::numeric_limits::epsilon() cannot be
accessed in Cuda. So I stick with the original implementation for now.
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313780388
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313780388
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779702
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779631
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779245
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
Review comment:
Comment added
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779293
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
Review comment:
Comment added
This is an automated message from the Apache Git Service.
To respond to the message, pleas
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779520
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779149
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
Review comment:
Fixed
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313778996
##
File path: src/operator/linalg_impl.h
##
@@ -1234,6 +1234,137 @@ LINALG_GPU_SYEVD_WORKSPACE_QUERY(DnDsyevd, double)
#endif // __CUDACC__
+ GESVD
+
+// CPU/GPU-versions of LAPACK function "gesvd"
+
+template inline
+void check_gesvd(const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V) {
+ // Any checking that helps user debug potential problems.
+ CHECK_LE(V.size(0), V.size(1))
+<< "The second to last dimension of A must be less or equal to the "
+<< "last dimension";
+ CHECK_EQ(UT.size(0), UT.size(1))
+<< "UT must be square matrix";
+ CHECK_EQ(V.size(0), L.size(0))
+<< "V, L have incompatible sizes";
+ CHECK_EQ(V.size(0), UT.size(0))
+<< "V, UT must have compatible sizes";
+}
+
+#define LINALG_CPU_GESVD(fname, DType) \
+template<> inline \
+void linalg_gesvd(const Tensor& UT, \
+ const Tensor& L, \
+ const Tensor& V, \
+ const Tensor& work, \
+ Stream *s) { \
+ check_gesvd(UT, L, V); \
+ DType lwork(0); \
+ MXNET_LAPACK_##fname(MXNET_LAPACK_ROW_MAJOR, V.size(0), V.size(1), \
Review comment:
Yes, you are right. I have removed the query.
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313779074
##
File path: src/operator/linalg_impl.h
##
@@ -1234,6 +1234,137 @@ LINALG_GPU_SYEVD_WORKSPACE_QUERY(DnDsyevd, double)
#endif // __CUDACC__
+ GESVD
+
+// CPU/GPU-versions of LAPACK function "gesvd"
+
+template inline
+void check_gesvd(const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V) {
+ // Any checking that helps user debug potential problems.
+ CHECK_LE(V.size(0), V.size(1))
+<< "The second to last dimension of A must be less or equal to the "
+<< "last dimension";
+ CHECK_EQ(UT.size(0), UT.size(1))
+<< "UT must be square matrix";
+ CHECK_EQ(V.size(0), L.size(0))
+<< "V, L have incompatible sizes";
+ CHECK_EQ(V.size(0), UT.size(0))
+<< "V, UT must have compatible sizes";
+}
+
+#define LINALG_CPU_GESVD(fname, DType) \
+template<> inline \
+void linalg_gesvd(const Tensor& UT, \
+ const Tensor& L, \
+ const Tensor& V, \
+ const Tensor& work, \
+ Stream *s) { \
+ check_gesvd(UT, L, V); \
+ DType lwork(0); \
+ MXNET_LAPACK_##fname(MXNET_LAPACK_ROW_MAJOR, V.size(0), V.size(1), \
+ UT.dptr_, UT.stride_, L.dptr_, V.dptr_, V.stride_, \
+ &lwork, -1); \
+ int ret(MXNET_LAPACK_##fname(MXNET_LAPACK_ROW_MAJOR, V.size(0), V.size(1), \
+ UT.dptr_, UT.stride_, L.dptr_, V.dptr_,
V.stride_, \
+ work.dptr_, static_cast(lwork))); \
+ CHECK_EQ(ret, 0) << #fname << " failed in lapack on cpu."; \
+}
+
+LINALG_CPU_GESVD(sgesvd, float)
+LINALG_CPU_GESVD(dgesvd, double)
+
+// Mangle temp storage requirements for DType and int into a single
Review comment:
Removed
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313778504
##
File path: src/operator/c_lapack_api.h
##
@@ -361,6 +382,26 @@ inline void flip(int m, int n, DType *b, int ldb, DType
*a, int lda) {
MXNET_LAPACK_CWRAP_SYEVD(ssyevd, float)
MXNET_LAPACK_CWRAP_SYEVD(dsyevd, double)
+ #define MXNET_LAPACK_CWRAP_GESVD(func, dtype) \
Review comment:
Added
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313778597 ## File path: src/operator/linalg.h ## @@ -191,6 +191,28 @@ int linalg_syevd_workspace_query(const Tensor& A, const Tensor& L, Stream *s = 0); + GESVD + +// CPU/GPU-versions of LAPACK function "gesvd". Please refer to the +// LAPACK documentation for further details. +// Note: +// - V is input and output parameter (overwritten by A) Review comment: Fixed This is an automated message from the Apache Git Service. To respond to the message, please log on to GitHub and use the URL above to go to the specific comment. For queries about this service, please contact Infrastructure at: us...@infra.apache.org With regards, Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313778173
##
File path: src/operator/c_lapack_api.h
##
@@ -242,6 +249,20 @@ inline void flip(int m, int n, DType *b, int ldb, DType
*a, int lda) {
#define MXNET_LAPACK_sgetrf LAPACKE_sgetrf
#define MXNET_LAPACK_dgetrf LAPACKE_dgetrf
+ #define MXNET_LAPACK_CWRAP_GESVD(prefix, dtype) \
+ inline int MXNET_LAPACK_##prefix##gesvd(int matrix_layout, int m, int n,
dtype* ut, \
Review comment:
The LAPACK_gesvd function interface differs in signature from the
MXNET_LAPACK-signature and have to be wrapped (as is stated
[here](https://github.com/apache/incubator-mxnet/blob/67191c4df9ba363605c59b332c1a8096573309e0/src/operator/c_lapack_api.h#L201)).
So this is basically a wrapper of LAPACK_gesvd.
I added some comments about how to use the LAPACK_gesvd. Its official
document can be found
[here](https://software.intel.com/en-us/mkl-developer-reference-c-gesvd).
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313710609
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
Review comment:
imho, X (square) overwritten by X L
This is an automated message from the Apache Git Service.
To respond
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313413734
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313412351
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313399704
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
Review comment:
Not yet. The public technical report (https://arxiv.org/pdf/1710.08717.pdf)
does not include details about svd.
This is an automated message from the Apache Git Service.
To respond to the message, please log on to GitHub and use the
URL above to go to the specific comment.
For queries about this service, please contact Infrastructure at:
us...@infra.apache.org
With regards,
Apache Git Services
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313397296
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313386657
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313386657
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313386657
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
[GitHub] [incubator-mxnet] hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
hzfan commented on a change in pull request #15795: [Numpy] Differentiable svd
URL: https://github.com/apache/incubator-mxnet/pull/15795#discussion_r313377236
##
File path: src/operator/numpy/linalg/np_gesvd-inl.h
##
@@ -0,0 +1,298 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one
+ * or more contributor license agreements. See the NOTICE file
+ * distributed with this work for additional information
+ * regarding copyright ownership. The ASF licenses this file
+ * to you under the Apache License, Version 2.0 (the
+ * "License"); you may not use this file except in compliance
+ * with the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing,
+ * software distributed under the License is distributed on an
+ * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+ * KIND, either express or implied. See the License for the
+ * specific language governing permissions and limitations
+ * under the License.
+ */
+
+/*!
+ * Copyright (c) 2017 by Contributors
+ * \file np_gesvd-inl.h
+ * \brief Function definition of the SVD Operator.
+ */
+#ifndef MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+#define MXNET_OPERATOR_NUMPY_LINALG_NP_GESVD_INL_H_
+
+#include
+#include
+#include "../../tensor/la_op.h"
+#include "../../tensor/la_op-inl.h"
+
+namespace mxnet {
+namespace op {
+
+struct GesvdVecSign {
+ template
+ MSHADOW_XINLINE static void Map(int i, int m, int n, DType* UT,
+ DType* V, int ldut, int ldv) {
+DType* vrow(V + i * ldv);
+DType maxval(fabs(vrow[0])), vval(0.0);
+int maxind(0);
+for (int i = 1; i < n; ++i) {
+ vval = fabs(vrow[i]);
+ if (vval > maxval) {
+maxval = vval;
+maxind = i;
+ }
+}
+if (vrow[maxind] < 0) {
+ DType* utcol(UT + i % m + (i / m) * ldut * m);
+ for (int i = 0; i < n; ++i) {
+vrow[i] = -vrow[i];
+if (i < m) {
+ utcol[i * ldut] = -utcol[i * ldut];
+}
+ }
+}
+ }
+};
+
+// (UT, L, V) = gesvd(A) [singular value decomposition]
+// - V can overwrite A
+// - Needs workspace (both DType and int), size of which is determined by a
+// workspace query
+struct gesvd {
+ template
+ static void op(const Tensor& A,
+ const Tensor& UT,
+ const Tensor& L,
+ const Tensor& V,
+ const OpContext& ctx,
+ const nnvm::NodeAttrs& attrs) {
+Stream *s = ctx.get_stream();
+if (A.dptr_ != V.dptr_) Copy(V, A, s);
+// From here on, we work on V only
+// Reserve workspace (size determined by query)
+int lwork(linalg_gesvd_workspace_query(UT[0], L[0], V[0], s));
+Tensor work = ctx.requested[0]
+ .get_space_typed(Shape1(lwork), s);
+// Loop over items in batch
+for (index_t i = 0; i < UT.size(0); ++i) {
+ linalg_gesvd(UT[i], L[i], V[i], work, s);
+}
+// Set signs in a deterministic way
+using namespace mxnet_op;
+Kernel::Launch
+ (s, V.size(0) * V.size(1), V.size(1), V.size(2),
+ UT.dptr_, V.dptr_, UT.stride_, V.stride_);
+ }
+};
+
+// (A) => (UT, L, V)
+template
+void NumpyLaGesvdForward(const nnvm::NodeAttrs& attrs,
+ const OpContext& ctx,
+ const std::vector& inputs,
+ const std::vector& req,
+ const std::vector& outputs) {
+ using namespace mshadow;
+ CHECK_EQ(inputs.size(), 1);
+ CHECK_EQ(outputs.size(), 3);
+ if (inputs[0].shape_.Size() == 0) {
+return;
+ }
+ MSHADOW_SGL_DBL_TYPE_SWITCH(outputs[0].type_flag_, OType, {
+mshadow::Stream *s = ctx.get_stream();
+laop::op(inputs[0].FlatToKD(s),
+ outputs[0].FlatToKD(s),
+ outputs[1].FlatToKD(s),
+ outputs[2].FlatToKD(s), ctx, attrs);
+ });
+}
+
+// Helper for gesvd_backward. See technical report for details
+template
+DType gesvd_back_helper_eps(DType* X);
+
+template<>
+MSHADOW_XINLINE float gesvd_back_helper_eps(float* X) {
+ return 1e-30;
+}
+
+template<>
+MSHADOW_XINLINE double gesvd_back_helper_eps(double* X) {
+ return 1e-100;
+}
+
+struct GesvdBackHelper_dV {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* L, int ldl,
+ DType* dA, int ldda) {
+const int offl(k * ldl);
+const int offda(k * m * ldda);
+DType denom(0.0);
+const DType eps(gesvd_back_helper_eps(dA));
+for (int i = 0; i < m; ++i) {
+ denom = L[offl + i];
+ if (denom < eps) denom = eps;
+ for (int j = 0; j < n; ++j) {
+dA[offda + i * ldda + j] /= denom;
+ }
+}
+ }
+};
+
+struct GesvdBackHelper_G1 {
+ template
+ MSHADOW_XINLINE static void Map(int k, int m, int n, DType* X, int ldx,
+ DType* L, int ldl) {
+const int offl(k * ldl);
+const int offx
