How to improve performance in a function that repeatedly defines and multiplies matrices

[DNF]

@Uranium238

# F1 = sqrt(det(Si*Sf*ρ^2*Wapp^-1))
# F1 = sqrt(det(Si*Sf*ρ^2)/det(W))
F1 = sqrt(det(Si*Sf*  (Wapp\ρ^2)))

To nitpick, ρ^2 * Wapp^-1 is not the same as Wapp \ ρ^2, but should be translated to ρ^2 / Wapp. Inside the determinant, it doesn’t matter much, though it leads to some floating point differences.

For performance (in addition to my doubts about the parallelization), there are many things that can be sped up, the most glaring one is perhaps

@Uranium238

for m in range(1,sysize)
        for mp in range(1,sysize)
            @. X11 = 0.0
            @. X12 = 0.0
            @. X21 = 0.0
            @. X22 = 0.0
            X11[m,:]=(J'*Si*J)[mp,:]*(-Bf[m,m])
            X12[m,:]=(J'*Bi*J)[mp,:]*(Bf[m,m])
            X21[m,:]=(J'*Si*J)[mp,:]*(Sf[m,m])
            X22[m,:]=(J'*Bi*J)[mp,:]*(-Sf[m,m])
            X = [X11 X12
                X21 X22]
            @. Y1 = 0.0
            @. Y2 = 0.0
            Y1[m,:] = Bf[m,m]*(J'*U*D)[mp,:] 
            Y2[m,:] = -Sf[m,m]*(J'*U*D)[mp,:]
            Y = [Y1
                Y2]
            g = im*tr(W\X)+    (transpose(W\V)*X*(W\V))[1]      -    (transpose(Y)*(W\V))[1]
            s = s+g[1]*gfcpart*T[m,mp]
        end
    end

where the products J' * Si * J etc. are performed over and over, instead of doing them once outside the loop. Also slicing along rows is slower than along columns, and the slices also allocate unnecessary temporary arrays, and the entire matrices X11, etc. are zeroed out for each iteration, even though only one row has been touched.

So there’s good news: at least GHT should be possible to speed up a lot