Граница множества Мандельброта

Тема в разделе "OpenCL", создана пользователем _qwe8013, 19 дек 2016.

  1. _qwe8013

    _qwe8013 Active Member

    Публикаций:
    2
    Регистрация:
    30 ноя 2016
    Сообщения:
    123
    Я смотрю, это будет первая тема в данном разделе.
    Понадобилось мне тут найти приближённо границу множества Мандельброта для [math]z^4+c[/math]. Создал шейдер:
    Код (C):
    1. float2 Mult(float2 z1,float2 z2)
    2. {
    3.     return (float2)(z1.x*z2.x-z1.y*z2.y, z1.y*z2.x+z1.x*z2.y);
    4. }
    5.  
    6. float2 DoIter(float2 z,float2 c)
    7. {
    8.     return Mult(Mult(Mult(Mult(z,z),z),z),z)+c;
    9. }
    10.  
    11. __kernel void CalcMandel(float Border,float Threshold,uint NumOfIter,uint Size,__global uchar *Data/*,__global uint *NumOfPoints*/)
    12. {
    13.     size_t _X=get_global_id(0);
    14.     size_t _Y=get_global_id(1);
    15.     if((_X<Size) && (_Y<Size))
    16.     {
    17.         float2 c=(float2)(Border*2.0f*((float)_X/(float)(Size-1)-0.5f),Border*2.0f*((float)(Size-_Y-1)/(float)(Size-1)-0.5f));
    18.         float2 z=(float2)(0,0);
    19.         for(uint i=0;(i<NumOfIter) && (length(z)<Threshold);i++)
    20.         {
    21.             z=DoIter(z,c);
    22.         }
    23.         uchar Val;
    24.         if(length(z)<Threshold)
    25.         {
    26.             Val=1;
    27.         }
    28.         else
    29.         {
    30.             Val=0;
    31.         }
    32.         Data[_X+Size*_Y]=Val;
    33.         barrier(CLK_GLOBAL_MEM_FENCE);
    34.         if(Val==1)
    35.         {
    36.             if((Data[_X+Size*_Y-1]==1) && (Data[_X+Size*_Y+1]==1))
    37.             {
    38.                 if((Data[_X+Size*(_Y-1)-1]==1) && (Data[_X+Size*(_Y-1)]==1) && (Data[_X+Size*(_Y-1)+1]==1))
    39.                 {
    40.                     if((Data[_X+Size*(_Y+1)-1]==1) && (Data[_X+Size*(_Y+1)]==1) && (Data[_X+Size*(_Y+1)+1]==1))
    41.                     {
    42.                         Val=0;
    43.                     }
    44.                 }
    45.             }
    46.         }
    47.         barrier(CLK_GLOBAL_MEM_FENCE);
    48.         Data[_X+Size*_Y]=Val;
    49.     }
    50.     else
    51.     {
    52.         barrier(CLK_GLOBAL_MEM_FENCE);
    53.         barrier(CLK_GLOBAL_MEM_FENCE);
    54.     }
    55. }
    Передаю в kernel параметры:
    Border=1.5
    Threshold=3.0
    NumOfIter=500
    Size=1000
    Data - массив байт.
    Жду на выходе границу, а получаю:
    [​IMG]
    на этой картинке чёрные точки соответствуют значению 1, а белые - 0.
    Что я делаю не так?
     

    Вложения:

    • Mandel.png
      Mandel.png
      Размер файла:
      30,7 КБ
      Просмотров:
      1.318
  2. _qwe8013

    _qwe8013 Active Member

    Публикаций:
    2
    Регистрация:
    30 ноя 2016
    Сообщения:
    123
    Добавлю все исходники.
     

    Вложения:

  3. _edge

    _edge Well-Known Member

    Публикаций:
    1
    Регистрация:
    29 окт 2004
    Сообщения:
    631
    Адрес:
    Russia
    Просьба также выкладывать EXE. Спасибо.
     
  4. _qwe8013

    _qwe8013 Active Member

    Публикаций:
    2
    Регистрация:
    30 ноя 2016
    Сообщения:
    123
    Вот EXE, сохраняет картинку в C:\temp\Mandel.bmp
     

    Вложения:

    • Mandel.zip
      Размер файла:
      12,7 КБ
      Просмотров:
      460
  5. murder

    murder Member

    Публикаций:
    0
    Регистрация:
    3 июн 2007
    Сообщения:
    628
    _qwe8013
    Хоть бы написал, что bmp сохраняется в C:\temp. У меня в системе например нет этой папки, пришлось создать. А вообще для сохранения в temp нужно использовать переменную среды %TEMP%.
     
  6. murder

    murder Member

    Публикаций:
    0
    Регистрация:
    3 июн 2007
    Сообщения:
    628
    Не совсем то но прикольно
    Код (C):
    1. float2 Mult(float2 z1,float2 z2)
    2. {
    3.    return (float2)(z1.x*z2.x-z1.y*z2.y, z1.y*z2.x+z1.x*z2.y);
    4. }
    5.  
    6. float2 DoIter(float2 z,float2 c)
    7. {
    8.    return Mult(Mult(Mult(Mult(z,z),z),z),z)+c;
    9. }
    10.  
    11. __kernel void CalcMandel(float Border,float Threshold,uint NumOfIter,uint Size,__global uchar *Data/*,__global uint *NumOfPoints*/)
    12. {
    13.    size_t _X=get_global_id(0);
    14.    size_t _Y=get_global_id(1);
    15.    if((_X<Size) && (_Y<Size))
    16.    {
    17.      float2 c=(float2)(Border*2.0f*((float)_X/(float)(Size-1)-0.5f),Border*2.0f*((float)(Size-_Y-1)/(float)(Size-1)-0.5f));
    18.      float2 z=(float2)(0,0);
    19.      for(uint i=0;(i<NumOfIter) && (length(z)<Threshold);i++)
    20.      {
    21.        z=DoIter(z,c);
    22.      }
    23.      uchar Val;
    24.      if(round(length(z)*0.2)==round(Threshold*0.2))
    25.      {
    26.        Val=1;
    27.      }
    28.      else
    29.      {
    30.        Val=0;
    31.      }
    32.      Data[_X+Size*_Y]=Val;
    33.      barrier(CLK_GLOBAL_MEM_FENCE);
    34.      barrier(CLK_GLOBAL_MEM_FENCE);
    35.    }
    36.    else
    37.    {
    38.      barrier(CLK_GLOBAL_MEM_FENCE);
    39.      barrier(CLK_GLOBAL_MEM_FENCE);
    40.    }
    41. }
    Mandel.png
     
  7. murder

    murder Member

    Публикаций:
    0
    Регистрация:
    3 июн 2007
    Сообщения:
    628
    Мне кажется, что это условие нужно запускать после того, как массив Data будет заполнен и результат выводить в другой массив (например Data2)
    Код (C):
    1.  
    2. barrier(CLK_GLOBAL_MEM_FENCE);
    3.         if(Val==1)
    4.         {
    5.             if((Data[_X+Size*_Y-1]==1) && (Data[_X+Size*_Y+1]==1))
    6.             {
    7.                 if((Data[_X+Size*(_Y-1)-1]==1) && (Data[_X+Size*(_Y-1)]==1) && (Data[_X+Size*(_Y-1)+1]==1))
    8.                 {
    9.                     if((Data[_X+Size*(_Y+1)-1]==1) && (Data[_X+Size*(_Y+1)]==1) && (Data[_X+Size*(_Y+1)+1]==1))
    10.                     {
    11.                         Val=0;
    12.                     }
    13.                 }
    14.             }
    15.         }
    16.         barrier(CLK_GLOBAL_MEM_FENCE);
    17.         Data2[_X+Size*_Y]=Val;
    Не?
     
    Mikl___ нравится это.
  8. _qwe8013

    _qwe8013 Active Member

    Публикаций:
    2
    Регистрация:
    30 ноя 2016
    Сообщения:
    123
    В #4 написано.
    Зачем Data2, если Data весь заполнен (для этого я и поставил barrier(CLK_GLOBAL_MEM_FENCE)).
     
  9. murder

    murder Member

    Публикаций:
    0
    Регистрация:
    3 июн 2007
    Сообщения:
    628
    Извиняюсь - прочитал немного про OpenCL теперь стало понятнее.

    Если втупую расчитывать каждую соседнюю точку - результат нормальный
    Код (C):
    1. float2 Mult(float2 z1,float2 z2)
    2. {
    3.    return (float2)(z1.x*z2.x-z1.y*z2.y, z1.y*z2.x+z1.x*z2.y);
    4. }
    5.  
    6. constant float Threshold;
    7. constant uint  Size;
    8. constant float Border;
    9. constant uint  NumOfIter;
    10.  
    11. uchar Mandel(size_t _X, size_t _Y){
    12.    float2 c=(float2)(Border*2.0f*((float)_X/(float)(Size-1)-0.5f),Border*2.0f*((float)(Size-_Y-1)/(float)(Size-1)-0.5f));
    13.    float2 z=(float2)(0,0);
    14.    for(uint i=0;(i<NumOfIter) && (length(z)<Threshold);i++)
    15.      z=Mult(Mult(Mult(Mult(z,z),z),z),z)+c;
    16.    if(length(z)<Threshold)
    17.      return 1;
    18.    else
    19.      return 0;
    20. }
    21.  
    22.  
    23. __kernel void CalcMandel(float _Border,float _Threshold,uint _NumOfIter,uint _Size,__global uchar *Data/*,__global uint *NumOfPoints*/)
    24. {
    25.    Threshold=_Threshold;
    26.    Size  =_Size;
    27.    Border  =_Border;
    28.   NumOfIter=_NumOfIter;
    29.    size_t _X=get_global_id(0);
    30.    size_t _Y=get_global_id(1);
    31.    Data[_X+_Size*_Y]=Mandel(_X,_Y);
    32.    if((Mandel(_X-1,_Y)==1)&&(Mandel(_X+1,_Y)==1))
    33.      if((Mandel(_X-1,_Y-1)==1)&&(Mandel(_X,_Y-1)==1)&&(Mandel(_X+1,_Y-1)==1))
    34.        if((Mandel(_X-1,_Y+1)==1)&&(Mandel(_X,_Y+1)==1)&&(Mandel(_X+1,_Y+1)==1))
    35.          Data[_X+_Size*_Y]=0;
    36. }
    1.png
    Я вот думаю val - это ведь private memory. И при синхронизации потоков в GPU просто не хватает памяти, чтобы сохранить значение переменной для каждого потока. Если сохранять промежуточные вычисления в __global массиве, то всё должно работать.
     
  10. _qwe8013

    _qwe8013 Active Member

    Публикаций:
    2
    Регистрация:
    30 ноя 2016
    Сообщения:
    123
    Так не должно быть, это - глюк.
    В вашем варианте производится больше вычислений.