X-ray scattering intensity for polymer simulation

Hello!

I have xyz data for a polymer that I simulated in lammps. I need to generate the x-ray scattering intensity profile for my system. Does anyone know how to do this? Here is a matlab script that ChatGPT kindly provided, but I don’t know if it’s right. I’m running another simulation that will allow me to compare with experimental data.

% Define atomic positions
positions = positions4; %(this is the dump file from lammps, brought into Matlab as a table and then table2array for just the xyz columns)
lambda = 0.15;  %X-ray wavelength in nm
L = 50; % Largest dimension of the polymer system in nm
theta_max = 30;  %Maximum half-scattering angle in degrees

% Define scattering vectors
q_min = 2 * pi / L;
q_max = (4 * pi / lambda) * sind(theta_max);
num_points = 10000; % number of q points
q_values = linspace(q_min, q_max, num_points);  % 1D array of q values

% Compute the structure factor

S_q = structure_factor(positions, q_values);

% Calculate the intensity
I_q = abs(S_q).^2;

% Convert q to 2*theta
two_theta = 4 * asind((q_values * lambda) / (4 * pi));

% Plot the intensity vs q
figure;
plot(two_theta, I_q);
xlabel('2theta');
ylabel('Intensity');
title('X-ray Scattering Intensity');
grid on;


% Function to compute the structure factor

function S_q = structure_factor(positions, q_values)
    N = size(positions, 1);
    S_q = zeros(size(q_values), 'like', 1 + 1i); % Initialize complex array

    for q_index = 1:length(q_values)
        q = q_values(q_index);
        sum_real = 0;
        sum_imag = 0;
        for j = 1:N
            pos = positions(j, :);
            qr = dot([q, 0, 0], pos); % q·r_j
            sum_real = sum_real + cos(qr);
            sum_imag = sum_imag + sin(qr);
        end
        S_q(q_index) = sum_real + 1i * sum_imag;
    end
end

For reference, what I need is found in this paper, figure 9: (DOI 10.1021/ma9810565)

A sample dump file is attached.

Thanks for your time.

AAPoly366dump.20500000.dump (429.3 KB)
Kind regards,
Sean

MAICoS has a module for small angle x-ray, it can return the structure factor or/and the intensity, see Small-angle X-ray scattering -

Simon

Isn’t the compute xrd just the right for that in LAMMPS?

I’ve never used it though.

Have a look at dynasor https://dynasor.materialsmodeling.org/

If you are interested in non-periodic structures you may use a code I have written some time ago GitHub - evoyiatzis/Jupyter-Notebooks: Jupyter notebooks for post-processing LAMMPS files
In the diffraction.py & diffraction_library.py files

Thank you, Simon!

Thanks for mentioning that, Germain! I didn’t even think to look at capabilities within lammps. I’ve got a simulation that just finished running where I tried it out.

Thanks! By non-periodic structures, do you mean something like non-periodic layer crystals (NPLs), in a polymer context? Because that is indeed what I am trying to characterize.

I have used this approach for computing the x ray pattern of metallic nanoparticles (see Eq. 5 in BJNANO - Atomistic insights into the morphological dynamics of gold and platinum nanoparticles: MD simulations in vacuum and aqueous media and references 80 & 81). Maybe you can have a look at the references and figure out if it is applicable to non periodic layer crystals.