Smart summation in FOR(mula) TRAN(slation)

Question

Well, today I'm here with this doubt...

I want to write this equation in Fortran:

Of course I can take the "classic" approach and write it like this:

 do i=1,N
        ac=0.0
        do j=i+1,M
            ac=ac+A(i,j)*B(j)
        enddo
        B(i)=(B(i)-ac)/A(i,i)
  enddo

But since I'm writing in Fortran, I would like to write it with an expression that looks like "more like the original" and, thus, compact. I was thinking in something like:

forall(i=1:N,j=i+1:M)
    B(i)=(B(i)- <MySummationExpression(A(i,j)*B(j))> )/A(i,i)
endforall

Expression that looks much more like the original. But the thruth is that I'm having a hard time finding a way to write the summation expression in a simple and compact way. Of course that I can write a function "real function summation(<expression>,<lower bound>, <upper bound>)", but since we're talking about Fortran I thought that there should be a simple (maybe intrinsic(?)) way of writing it.

So there is a compact way to write that expression, or I have to take the uglier way (two explicit loops)?

EDIT: In the actual code x is a two dimensional array, with one solution in each column. So using the intrinsic function sum that, up to here, seems like a great idea (as @alexander-vogt shown in his answer) leads to almost the same "compacticity" of the code:

do j=1,size(B,2)
        do i=nA,1,-1
            B(i,j)=(B(i,j)-sum(A(i,i+1:nA)*B(i+1:nA,j)))/A(i,i)
        enddo
 enddo

Answer 1

What about:

do i=1,N
  B(i) = (B(i) - sum( A(i+1:M,i) * B(i+1:M) )) / A(i,i)
enddo  

(Please notice that I changed the indexing of matrix A, Fortran is column major! )

As we are re-using B for the result, this loop has to be done in ascending order. I'm not sure it's possible to do this using forall (isn't it up to the compiler to choose how to proceed? - see Fortran forall restrictions).

Writing the result into a new vector C does not overwrite B and can be executed in any order:

forall ( i=1:N )
  C(i) = (B(i) - sum( A(i+1:M,i) * B(i+1:M) )) / A(i,i)
endforall