خوارزمية أصدقاء الأصدقاء المكتوبة بلغة بايثون يجب أن تكون في فورتران 90/95

Question

أحاول كتابة الكود الخاص بي لخوارزمية "أصدقاء الأصدقاء".تعمل هذه الخوارزمية على مجموعة من نقاط البيانات ثلاثية الأبعاد وترجع عدد "الهالات" في مجموعة البيانات.كل هالة عبارة عن مجموعة من النقاط التي تكون المسافة بينها أقل من طول الارتباط، b، المعلمة الوحيدة للبرنامج.

الوصف الخوارزمي:تحتوي خوارزمية FOF على معلمة واحدة مجانية تسمى طول الارتباط.يُطلق على أي جسيمين تفصل بينهما مسافة أقل من أو تساوي طول الارتباط اسم "الأصدقاء".يتم بعد ذلك تعريف مجموعة FOF من خلال مجموعة الجسيمات التي يرتبط فيها كل جسيم داخل المجموعة بكل جسيم آخر في المجموعة من خلال شبكة من الأصدقاء.

اضبط عداد مجموعة FOF j=1.

• لكل جسيم، n، غير مرتبط بعد بأي مجموعة:

• قم بتعيين n للمجموعة j، وقم بتهيئة قائمة أعضاء جديدة، mlist، للمجموعة j مع الجسيم n كمدخل أول،

• بشكل متكرر، لكل جسيم جديد p في القائمة:

• ابحث عن جيران p الذين يقعون على مسافة أقل من أو تساوي طول الارتباط، وأضف إلى القائمة تلك التي لم يتم تعيينها بالفعل للمجموعة j،
• سجل قائمة ml للمجموعة j، اضبط j=j+1.

هذه هي محاولتي لترميز الخوارزمية.اللغة الوحيدة التي أشعر بالارتياح للقيام بها هي لغة بايثون.ومع ذلك، أحتاج إلى كتابة هذا الرمز بلغة فورتران أو جعله أسرع.آمل حقًا أن يساعدني شخص ما.

أولاً أقوم بإنشاء مجموعة من النقاط التي يجب أن تحاكي وجود 3 هالات:

import random
from random import *
import math
from math import *
import numpy
from numpy import *
import time

points = 1000

halos=[0,100.,150.]

x=[]
y=[]
z=[]
id=[]
for i in arange(0,points,1):
   x.append(halos[0]+random())
   y.append(halos[0]+random())
   z.append(halos[0]+random())
   id.append(i)

for i in arange(points,points*2,1):
   x.append(halos[1]+random())
   y.append(halos[1]+random())
   z.append(halos[1]+random())
   id.append(i)

for i in arange(points*2,points*3,1):
   x.append(halos[2]+random())
   y.append(halos[2]+random())
   z.append(halos[2]+random())
   id.append(i)

ثم قمت بترميز خوارزمية FOF:

  x=array(x)
  y=array(y)
  z=array(z)
  id=array(id)

  t0 = time.time()                         

  id_grp=[]
  groups=zeros((len(x),1)).tolist()
  particles=id
  b=1 # linking length
  while len(particles)>0:
  index = particles[0]
  # remove the particle from the particles list
  particles.remove(index)
  groups[index]=[index]
  print "#N ", index
  dx=x-x[index]
  dy=y-y[index]
  dz=z-z[index]
  dr=sqrt(dx**2.+dy**2.+dz**2.)
  id_to_look = where(dr<b)[0].tolist()
  id_to_look.remove(index)
  nlist = id_to_look
  # remove all the neighbors from the particles list
  for i in nlist:
        if (i in particles):
           particles.remove(i)
  print "--> neighbors", nlist
  groups[index]=groups[index]+nlist
  new_nlist = nlist
  while len(new_nlist)>0:
          index_n = new_nlist[0]
          new_nlist.remove(index_n)
          print "----> neigh", index_n
          dx=x-x[index_n]
          dy=y-y[index_n]
          dz=z-z[index_n]
          dr=sqrt(dx**2.+dy**2.+dz**2.)
          id_to_look = where(dr<b)[0].tolist()
          id_to_look = list(set(id_to_look) & set(particles))
          nlist = id_to_look
          if (len(nlist)==0):
             print "No new neighbors found"
          else:
             groups[index]=groups[index]+nlist
             new_nlist=new_nlist+nlist
             print "------> neigh-neigh", new_nlist
             for k in nlist:
               particles.remove(k)

في النهاية ينتهي الأمر بقائمة من الهالات الموجودة في القائمة groups

هذا الجزء من الكود خارج الموضوع قليلاً ولكن أعتقد أنه سيكون من الجيد أن أعرضه لك.أقوم بشكل أساسي بحذف جميع المجموعات التي لا تحتوي على جزيئات، وفرزها وفقًا لعدد الجزيئات وإظهار بعض الخصائص.

  def select(test,list):
  selected = []
  for item in list:
    if test(item) == True:
      selected.append(item)
  return selected

  groups=select(lambda x: sum(x)>0.,groups)
  # sorting groups
  groups.sort(lambda x,y: cmp(len(x),len(y)))
  groups.reverse()

  print time.time() - t0, "seconds"

  mass=x
  for i in arange(0,len(groups),1):
    total_mass=sum([mass[j] for j in groups[i]])
    x_cm = sum([mass[j]*x[j] for j in groups[i]])/total_mass
    y_cm = sum([mass[j]*y[j] for j in groups[i]])/total_mass
    z_cm = sum([mass[j]*z[j] for j in groups[i]])/total_mass
    dummy_x_cm = [x[j]-x_cm for j in groups[i]]
    dummy_y_cm = [y[j]-y_cm for j in groups[i]]
    dummy_z_cm = [z[j]-z_cm for j in groups[i]]
    dummy_x_cm = array(dummy_x_cm)
    dummy_y_cm = array(dummy_y_cm)
    dummy_z_cm = array(dummy_z_cm)
    dr = max(sqrt(dummy_x_cm**2.+dummy_y_cm**2.+dummy_z_cm**2.))
    dummy_x_cm = max(dummy_x_cm)
    dummy_y_cm = max(dummy_y_cm)
    dummy_z_cm = max(dummy_z_cm)
    print i, len(groups[i]), x_cm, y_cm, z_cm,dummy_x_cm,dummy_y_cm,dummy_z_cm

Answer 1

أعتقد أنه من غير المستحسن أن تبدأ في تعلم فورتران على أمل أن يكون الكود الناتج أسرع من تطبيقك الحالي.قد يكون الأمر كذلك في النهاية، ولكن أعتقد أنه من الأفضل أن تقوم بتنفيذ Python الخاص بك بأسرع ما يمكن قبل التفكير في التنفيذ بلغة أخرى، خاصة لغة أجنبية.

أنا أكتب Fortran، وأعتقد شخصيًا أن أداءه مزعج في جميع أنحاء Python، لكن الأشخاص الذين يعرفون هذه الأشياء يقدمون حججًا مقنعة مفادها أن Python + SciPy + Numpy يمكن، إذا تم تصميمها بعناية، أن تضاهي Fortran من حيث السرعة في النواة الحسابية للعديد من العلوم/الهندسية. البرامج.لا تنس أنك لم تقم بتحسين لغة Python الخاصة بك حتى يتم تشغيل جميع النوى الموجودة على جهاز الكمبيوتر الخاص بك باللون الأحمر.

لذا:

الأول - الحصول على تطبيق عملي في بايثون.

2 - اجعل التنفيذ في أسرع وقت ممكن.

إذا (أحرف كبيرة لأنها كبيرة "إذا") لا يزال الرمز غير سريع بما فيه الكفاية وتكون تكلفة/فائدة ترجمته إلى لغة مجمعة مناسبة، ففكر في اللغة المترجمة التي تريد ترجمتها إليها.إذا كنت تعمل في مجال يستخدم فيه لغة Fortran على نطاق واسع، فتعلم لغة Fortran بكل الوسائل، ولكنها لغة متخصصة وقد يفيد مسيرتك المهنية أكثر أن تتعلم لغة C++ أو أحد أقاربها.

يحرر (طويل جدًا بحيث لا يمكن وضعه في مربع التعليق)

لماذا تضللنا في سؤالك؟لقد ذكرت أن اللغة الوحيدة التي تشعر بالارتياح معها هي لغة بايثون، والآن تقول أنك تعرف لغة فورتران.أعتقد أنك يجب أن تكون غير مريح معها.ومن تعليقك، يبدو أن المساعدة التي تحتاجها حقًا هي جعل تنفيذ Python أسرع؛ عرض جانبي بوب وقد عرضت بعض النصائح.خذ ذلك بعين الاعتبار، ثم قم بالتوازي.

Answer 2

A pointer towards a more efficient algorithm. If I am not mistaken you are comparing a point to every other point to see if any are closer than the linking length. For large numbers of points there are quicker ways of finding a near neighbours - spatial indexing and KD trees off the top of my head but doubtless there are other methods too that will work for you.

Answer 3

If you have a modern graphics card you can paralelize to hundreds of processors (depending on your card) in Python code using PyOpenCL.

You can investigate to see if the algoritm FoF is implemented inside this voidfinder F90 code

You can define distance as the squared distance to avoid the use of sqrt() and use x*x instead of x**2...