Fft python - analisi periodicità 2d gruppo atomi

[gaspa871]

Ciao. Mi sono iscritto al forum perchè ho bisogno di aiuto.
Sicuramente il mio problema per voi è banale e in realtà per questo vi scrivo.

Attualmente sono in tesi alla magistrale in scienza dei materiali a padova. Nella mia tesi mi occupo di simulazioni di dinamica molecolare. Una delle analisi che devo fare riguarda capire l'ordine di un gruppo di atomi di zolfo adsorbiti su una superficie ideale di oro 111. Devo verificare che mantengano o meno nel corso della simulazione un ordine esagonale. Ad ogni istante conosco le coordinate xyz degli zolfi adsorbiti.
Quello che voglio fare è quello che normalmente da sempre si fa in fisica per indigare la struttura della materia, cioè una trasformata di fourrier che a partire dal reticolo diretto definisce un preciso reticolo reciproco. In pratica voglio calcolare il fattore di struttura o l'intensità dei raggi difratti. Una cosa simile a quella che voglio fare è riportata in questo articolo :

http://www.jim.or.jp/journal/e/pdf3/51/04/664.pdf

Mi sono scritto in python uno scriptino che fa questa analisi, ma è molto lento e preferirei gli algoritmi fft implementati in numpy. Il problema è che non riesco a capire quale algoritmo utilizzare e come utilizzarli. I miei input sono le coordinate atomiche e vorrei avere come output il fattore di struttura o l'intensità dei raggi difratti. Non capisco come devo passare l'input a queste procedura per fare la trosfarmata nel modo riportato nell'articolo sopra citato.
Mi scuso per il disturbo e vi ringrazio per l'aiuto.

Buona serata,
Piero.

[xdom="Seneca"]Sposto la discussione in Matematica per l'Economia e per le Scienze Naturali.[/xdom]

[yoshiharu]

"gaspa87":
Mi sono scritto in python uno scriptino che fa questa analisi, ma è molto lento e preferirei gli algoritmi fft implementati in numpy. Il problema è che non riesco a capire quale algoritmo utilizzare e come utilizzarli. I miei input sono le coordinate atomiche e vorrei avere come output il fattore di struttura o l'intensità dei raggi difratti. Non capisco come devo passare l'input a queste procedura per fare la trosfarmata nel modo riportato nell'articolo sopra citato.

Ciao,
se non ho capito male, il fattore di struttura e'
[tex]\sum_{p\in P} e^{i \underline k\cdot \underline x_p}[/tex]

laddove $P$ e' l'insieme dei nodi del reticolo di atomi.
In questo caso, piu' che parlare di trasformata, io imposterei la cosa in maniera un pochino diversa.
Tanto per essere definiti, qualcosa del tipo (mi sembra di aver capito che sei nel caso 2D):

from numpy import *

def mapExp(a,k):
    return exp(1j * k * asarray(a))

def structFactor(a,b,k):
    return dot(mapExp(a,k[0]),mapExp(b,k[1]))

e il fattore di struttura sarebbe

structFactor(lista_ascisse,lista_ordinate,[k_x,k_y])

Praticamente ottieni la somma come prodotto scalare dei due vettori
[tex]\left(e^{i k_x x_0},\dots,e^{i k_x x_{N-1}}\right) \qquad \left(e^{i k_y y_0},\dots,e^{i k_y y_{N-1}}\right)[/tex]
Il codice assume che tu abbia le coordinate in liste (o altre iterabili) python.
Puoi generalizzare in piu' dimensioni ricordandoti che il prodotto tra due array di numpy viene fatto elemento per elemento.

In alternativa puoi sempre scriverti una estensione (vale a dire un modulo scritto in C), mettendo li' tutte le cose piu' "bisognose" dal punto di vista della CPU. Pero' effettivamente con numpy e' molto piu' semplice