# Finite differences with Julia

Dear Julia users,
Is there a way to improve performance of this kind of code without devectorizing (use loops). All variables are 2d arrays. I solve a PDE problem with finite differences method. Thanks for your help

``````996988748                divV .=  diff(Vx,dims=1)/dx
992096000                divV .+= diff(Vy,dims=2)/dy
-
992096000                Exxc .= diff(Vx,dims=1)/dx
495600000                Exxc .-= 1/2*divV
992096000                Eyyc .= diff(Vy,dims=2)/dy
495600000                Eyyc .-= 1/2*divV
-
507248000                Exyv  .= diff(Vx_exp,dims=2)
0                Exyv ./= dy
507248000                Exyv .+= diff(Vy_exp,dims=1)
0                Exyv ./= dx
0                Exyv .*= 0.5
-
1984192000         @views Exyc  .= 0.25 * ( Exyv[1:end-1,1:end-1]
-                                + Exyv[2:end  ,1:end-1]
-                                + Exyv[1:end-1,2:end  ]
-                                + Exyv[2:end  ,2:end  ])
``````

Don’t use `diff`, that creates temporary arrays. Don’t use multiple lines of `.+=`, `./=`, etcetera — you want a single fused loop for cache efficiency. Don’t use `*` and `+` and `/` when you want `.*` and `.+` and `./` in order to fuse the loops.

(But honestly, this is the sort of thing that is probably more readable if you just write a loop, and there is no performance reason not to do so in Julia.)

7 Likes

Thanks for your answer, i rewrite the code following your advices and it is better. With loops the code is faster. If i understood well, in julia the fortran-like syntax is better than the
matlab-like syntax. Unfortunately, the code i wanted to translate is in matlab I really like Julia as a language but fortran is still faster for what i am doing (for now).

``````  1024000         @views divV .= ((f.Vx[2:end,:] .- f.Vx[1:end-1,:])./dx
-                       .+(f.Vy[:,2:end] .- f.Vy[:,1:end-1])./dy)
-
1024000         @views Exyv .= 0.5 .*( (Vx_exp[:,2:end] .- Vx_exp[:,1:end-1])./dy
-                          .+  (Vy_exp[2:end,:] .- Vy_exp[1:end-1,:])./dx)
-
1024000         @views Exyc  .= 0.25 .* (  Exyv[1:end-1,1:end-1]
-                                 .+ Exyv[2:end  ,1:end-1]
-                                 .+ Exyv[1:end-1,2:end  ]
-                                 .+ Exyv[2:end  ,2:end  ])
-
0                Eii2  .= 0.5 .* (Exxc.^2 .+ Eyyc.^2) .+ Exyc.^2
``````

Check out this for more information. One example is on a reaction-diffusion equation.

2 Likes

Thank you very much for help, the last link will help
My code is here https://github.com/pnavaro/M2Dpt.jl
I will try to beat fortran and matlab during next few days.

I agree this may be clearer & faster with loops, or functions containing the loops. But sometimes it’s good to have something which matches the old code almost line-for-line.

You may want to try `@uviews` from UnsafeArrays.jl for creating lots of views you’ll only use once.