文件交换Pick of the Week

我们最好的用户提交

高效的面向对象的Kronecker产品操作

发布者Sean de Wolski那

31意见（过去30天）|Likes|1条评论
Sean 's pick this week is KronProd 经过 Matt J 。

Background

This one is old and very dear to my heart. Some time in the 2009 era, as I was starting grad school, I needed help calculating subpixel resolution displacement vectors between two CT image volumes. The displacement vector calculation was one of the main underlying pieces of most of my grad research. This was before imregcorr existed, though it still does not support 3D which was my requirement.
我会设法获得一个2D变体工作，但需要3D。Matt回复了我的新闻阅读帖子，讨论算法及其在3D中的扩展，并指出这比使用更有效 凯龙 directly. Performance was important as I'd have to compute over four billion displacement vectors over the next year.
Matt提供了许多算法，并重载了一堆运营商来与之合作 KronProd 类并使用自定义索引和大小配置它。2009年，我不明白如何定制索引 子反馈 / subsasgn./尺寸 重载工作，所以它看起来像魔法。凭借现在这些工作的知识，以及该时代的绩效影响，我可能会将核心算法拉到函数以便因为我的使用是无状态的，即一个新的 KronProd 每一次。
Here are a few simple operations with a KronProd 。We'll create one from a three arrays.
kp = KronProd({magic(3), [randi(10,[3 4]), zeros(3,1)], randn(9, 15)})
kp =
8.1×225 KronProd array with properties: opset: {[3×3 double] [3×5 double] [9×15 double]} opinds: [1 2 3] numops: 3 eyemask: [0 0 0] domainset: [3 5 15] maxdim: 3 scalarcoeff: 1 scalarcumprods: [1 1 1] scalarmask: [0 0 0] domainsizes: [3 5 15] rangesizes: [3 3 9]
现在我们可以使用它，就像你从中获得的完全扩张 凯龙 功能。
NZEROS = NNZ（KP）
纳泽斯= 14580.
Or a matrix operation where I get back KronProds.
[u, s, v] = svd(kp)
u =
81×81 KronProd阵列具有属性：{3×1 Cell} Opinds：[11] Numops：3 Eyemask：[0 0 0] Domainset：[3×1双] MaxDim：3 Scalarcoeff：1 ScalArcumProds：[1 1 1] Scalarmask：[0 0 0]域：[3 3 9]范围：[3 3 9]
s =
8.1×225 KronProd array with properties: opset: {3×1 cell} opinds: [1 2 3] numops: 3 eyemask: [0 0 0] domainset: [3×1 double] maxdim: 3 scalarcoeff: 1 scalarcumprods: [1 1 1] scalarmask: [0 0 0] domainsizes: [3 5 15] rangesizes: [3 3 9]
v =
225×225 KronProd阵列具有属性：Opset：{3×1 Cell} Opinds：[11] Numops：3 Eyemask：[0 0 0] Domaintet：[3×1双] MaxDim：3 Scalarcoeff：1 ScalArcumProds：[1 1 1] Scalarmask：[0 0 0]域：[3 5 15]范围：[3 5 15]
马特超负荷 子反馈 so when you index into the KronProd, you get back a KronProd or just matrix.
kkp = kp（1）
KKP =
3×3 KronProd阵列具有属性：Opset：{[3×3双]} Opinds：1 Numops：1 Eyemask：0 Domainset：3 MaxDim：1 ScalarCumProds：1 ScalarmumProds：1 Scalarmask：0域：3范围：3
km = kp {1}
km = 3×3
8.16.3.5 7 4 9 2
他没有过载 subsasgn. though, so we can't do the opposite.
尝试
kp{1} = rand(3)
catch我
fprintf(2, ME.message)
结尾
Unable to perform assignment because brace indexing is not supported for variables of this type.
他也过度了 尺寸 so the appearance of size is that of the full Kronecker product. Under the hood, however, the KronProd 实际上是一个标量，它只是 出现 成为扩展的大小。
sz =尺寸
SZ = 1×2
8.1225
Starting in MATLAB R2021b, we now have an easier and more formal means for controlling indexing and creating scalar objects using 模块化索引 。
Right now 子反馈 过载（），{}和点（。）。它超载了点，因为它必须用旧方法，所以我们不需要用更新的索引来做到这一点。所以类定义线需要：
classdefKronprod
Now, we'd need to implement the abstract methods for parenAssign, parenReference, braceAssign, braceReference 。然后，这些方法将在受保护的方法块中使用与今天使用的类似逻辑，沿着这些线路的东西：
功能B = braceReference(obj, indexOp)
inds = m.opinds（Indexop（1）.indices {1}）;
b = m.opset（Inds）;
如果isnum(inds)
b = b {：};
结尾
结尾

注释

马特，如果我们在现实生活中见面，我欠你啤酒！
试一试，让我们知道你的想法 here or leave a comment for Matt.
|
Category:
挑选
  • 打印
  • 发送电子邮件

也可以看看

注释

To leave a comment, please clickhereto sign in to your MathWorks Account or create a new one.