Hi Tim,
Many thanks for the suggestion, unfortunately it did not fix the issue; see
below updated code + timing/allocation.
All the best, Christoph
function evalSiteFunction!(sp::SiteLinear, geom::tbgeom,
Y::Array{Float64,2}, Frc::Array{Float64,2})
# extract dimension information (d1, d2 \in \{2, 3\})
d1, d2, nneig = size(sp.coeffs)
# allocate some variables
kX = 0::Int; jX = 0::Int
# assert type of the index set over which we are looping
J = geom["iMM"]::Array{Int, 1}
for jX in J
# initialise force to 0
Frc[:,jX] = 0.0
# loop over neighbouring sites
for n in 1:nneig
# get the X-index of the current neighbour
kX = geom.I[ geom.X2I[1,jX] + sp.stencil[1,n], geom.X2I[2,jX] +
sp.stencil[2,n] ]
# evaluate the term (devectorized for performance)
# f[:] += slice(sp.coeffs, :, :, n) * Y[:, kX]
for i = 1:d1, j = 1:d2
Frc[i,jX] += sp.coeffs[i, j, n] * Y[j, kX]
end
end
end
end
@time TB.evalSiteFunction!(frcMM, geom, Y, Frc)
elapsed time: 0.549490861 seconds (166123912 bytes allocated, 21.85% gc time)
1 abstractarray.jl; checkbounds; line: 62
1 array.jl; getindex; line: 247
420 task.jl; anonymous; line: 340
420 .../IJulia/src/IJulia.jl; eventloop; line: 123
420 ...rc/execute_request.jl; execute_request_0x535c5df2; line: 140
420 loading.jl; include_string; line: 97
420 profile.jl; anonymous; line: 14
45 ...qmmm/TBmultiscale.jl; evalSiteFunction!; line: 198
1 .../lib/julia/sys.dylib; +; (unknown line)
1 array.jl; getindex; line: 247
1 ...qmmm/TBmultiscale.jl; evalSiteFunction!; line: 202
372 ...qmmm/TBmultiscale.jl; evalSiteFunction!; line: 203
3 array.jl; getindex; line: 247
8 array.jl; setindex!; line: 308
