To calculate the structure factor O(N), where N is the number
of atoms, you can utilize S(\vec q)= | F(\vec q) |^2 where
the form factor amplitude is defined as
F(\vec q)=sum_i f_i exp(i \vec q\dot \vec r_i )
and f_i denote the scattering length of the i'th atom which
acounts for the details of the radiation-matter interaction.
Hence each process for each \vec q have to sum
f_i cos(q r_i) and f_i sin(q r_i) over all locally owned
atoms and store them in two arrays. Afterwards do
an MPI_ALLreduce(.., MPI_SUM , .. ) on the array. Then
the amplitude is known globally and cheap to square,
afterwards you can bin the structure factor after
\vec q| and average.
\vec q can only be points on the reciprocal lattice
corresponding to the box due to periodic boundary
conditions. For a fixed rectangular box (Lx,Ly,Lz):
\vec q (nx,ny,nz) = (nx 2Pi/Lx, ny 2Pi/Ly, nz 2Pi/Lz),
similar can be defined for fixed triclinic boxes.
Here FFTs is one option, but a fast and transparent
approach is to realise for a given atom and q vector \vec r_i \dot \vec q
is of the form nx a + ny b + nz c. Since cartesian directions factorise
you are left with the task of calculating cos(nx a) and sin(nx a) for
nx=0, 1, .., M. You can do nx=0,1 explicitly, and then express the case
=2 recursively with products and sums of previous values. In particlar,
for even nx cos(nx a) can be expressed with cos(a nx/2) and sin(a nx/2)
which are already known. For odd nx cos(nx a) can be expressed though
cos(a),sin(a), cos((nx-1)a) and sin((nx-1)a) which have already been
computed. That can be used to derive all the cos, sin terms up to
q><qmax evaluating only 6 trigonometric functions for each atom.
That should produce an algorithm which is O( N M^3) where N is number
of atoms and M^3 is the number of q vectors you want to evaluate.