你们好,我想将MATLAB代码改为Julia,但是最后得到结果不一致,其中偏差在矩阵运算那边一块。
我将某些运算提取出来单独运行,结果不一样,代码都贴在下面。
希望各位帮帮忙,感谢感谢。
Julia程序
A = [
[1.0 -1 2 -1]
[2 -2 3 -3]
[1 1 1 0]
[1 -1 4 3]
];
for k = 2:3
A[k+1:4,k] = A[k+1:4, k] ./ A[k, k]
A[k+1:4, k+1:4] = A[k+1:4, k+1:4] .- A[k+1:4, k] .* A[k,k+1:4]
end
A
Julia运行结果是:
4×4 Array{Float64,2}:
1.0 -1.0 2.0 -1.0
2.0 -2.0 3.0 -3.0
1.0 -0.5 2.5 1.5
1.0 0.5 2.2 1.2
MATLAB程序
A = [
[1.0 -1 2 -1]
[2 -2 3 -3]
[1 1 1 0]
[1 -1 4 3]
]
for k = 2:3
A(k+1:4,k) = A(k+1:4, k) ./ A(k, k);
A(k+1:4, k+1:4) = A(k+1:4, k+1:4) .- A(k+1:4, k) .* A(k,k+1:4);
end
A
MATLAB运行结果:
A =
1.00000 -1.00000 2.00000 -1.00000
2.00000 -2.00000 3.00000 -3.00000
1.00000 -0.50000 2.50000 -1.50000
1.00000 0.50000 1.00000 6.00000
楼主你好。关于这个计算式,
A[k+1:4, k] .* A[k,k+1:4]
matlab和julia的理解方式不同。
julia理解成两个一维向量逐像素相乘,仍得到一个一维向量。
而matlab理解成一个行向量和列向量相乘,得到了一个二维的矩阵。
你可以把它们打印出来,就能看到区别了。
1 个赞
原因见楼上。
把第二个向量转置一下就行了。
A = [
[1.0 -1 2 -1]
[2 -2 3 -3]
[1 1 1 0]
[1 -1 4 3]
];
for k = 2:3
A[k+1:4,k] = A[k+1:4, k] ./ A[k, k]
A[k+1:4, k+1:4] = A[k+1:4, k+1:4] .- A[k+1:4, k] .* A[k,k+1:4]' # <== 这里转置了
end
println("A = $A")
julia> A
4×4 Array{Float64,2}:
1.0 -1.0 2.0 -1.0
2.0 -2.0 3.0 -3.0
1.0 1.0 1.0 0.0
1.0 -1.0 4.0 3.0
julia> k=2
2
julia> A[k,k+1:4]
2-element Array{Float64,1}:
3.0
-3.0
julia> A[k+1:4, k] .* A[k,k+1:4]
2-element Array{Float64,1}:
3.0
3.0
julia> A[k,k+1:4]'
1×2 LinearAlgebra.Adjoint{Float64,Array{Float64,1}}:
3.0 -3.0
julia> A[k+1:4, k] .* A[k,k+1:4]'
2×2 Array{Float64,2}:
3.0 -3.0
-3.0 3.0
1 个赞
我来补充一下吧,在Julia中有一种操作叫做广播(broadcast),它可以把某个函数应用到容器里的每个元素上。而a.*b
是一种语法糖(简便写法的意思),相当于广播了*(a,b)
函数(在Julia里,二元运算符自身是一个二元函数):
broadcast(*,a,b)
对于Julia中的Array(多维数组,类似MATLAB中的矩阵)来说,广播有预定义的特殊含义。假设a和b都是一维的,且两个都是1*N或者都是N*1,那么广播的结果是形状和输入元素一样的数组,对应位置的元素是a和b对应位置进行*
运算的结果:
julia> a=[1,2,3]
3-element Array{Int64,1}:
1
2
3
julia> b=[4,5,6]
3-element Array{Int64,1}:
4
5
6
julia> broadcast(*,a,b)
3-element Array{Int64,1}:
4
10
18
当然,其中一个也可以是标量(或者0维数组),这时类似数乘向量的效果。
另一个需要注意的是,Julia中数组赋值时的形状也必须匹配。例如:
julia> A=[1 2;3 4]
2×2 Array{Int64,2}:
1 2
3 4
julia> A[:,1]=[5,7]
2-element Array{Int64,1}:
5
7
julia> A[:,1].=0
2-element view(::Array{Int64,2}, :, 1) with eltype Int64:
0
0
如果要将第一列刷成零,必须使用.=
运算符赋值。直接赋值会发生ArgumentError
。
4 个赞
razor
7
谢谢。那么向量化 和 广播 区别是什么呢?对这两词有点困惑。