Проблема с созданием Mandelbrot Задание графика в C

Question

Я пытаюсь создать набор Мандельброта в коде C; моим выходом будет файл данных, один столбец реальных частей и один столбец мнимых частей, которые будут нанесены на график в Арганде или комплексной плоскости. У меня есть весь сложный математический материал, определенный в заголовке complex.h, и я использую структуру complex с double R для реальной части и double I для мнимой. Я пытаюсь зацикливать значения для dz и обновлять z iter = 80 раз или пока точка не окажется за пределами определенного радиуса. dz мнимая, по существу, на диапазоне [(dzrmin + i * dzimin), (dzrmax + i * dzimax)]. z - текущее комплексное число, обновляемое как z^2 = dz + z. У меня есть функция csum(), чтобы суммировать два комплексных числа и функцию csquare(), чтобы правильно квадратировать комплексное число. Вот весь мой код

#include <stdio.h> 
    #include <stdlib.h> 
    #include <math.h> 
    #include "complex.h" 

    int main(void) 
    { 
     double dzrmin, dzrmax, dzimin, dzimax, dzr,dzi, r0,i0, maxrad; 
     int i,j,k,iter, Ndr,Ndi; 
     complex z, dz; 
     dzrmin = -2.1; 
     dzrmax = 0.6; 
     dzimin = -1.2; 
     dzimax = 1.2; 
     Ndr = 200 ; 
     Ndi = 180; 
     dzr = (dzrmax - dzrmin)/Ndr; 
     dzi = (dzimax - dzimin)/Ndi; 
     r0 = 0.0; 
     i0 = 0.0; 
     z.R = r0; 
     z.I = i0; 
     dz.R = dzrmin; 
     dz.I = dzimin; 
     maxrad = 2.0; 
     iter = 80 ; 
     for(i = 0; i < Ndr; i++) 
     { 
        for(j = 0; j < Ndi; j++) 
        { 

         for(k = 0; k< iter; k++) 
         { 
        printf("%.6lf %.6lf\n", z.R, z.I); 
        z = csum(csquare(z), dz) ; 

        if(cmag(z) > maxrad) k = iter ; 

       } 

       dz.I += dzi ; 
      } 
      dz.R += dzr; 
      z.R = r0; 
      z.I = i0; 
     } 
     return 0; 
    } 

и мой файл заголовка complex.h

#include <stdlib.h> 
#include <math.h> 
typedef struct complexnumber{ double R; double I ; } complex ; 
double cmag(complex z)                                    
{ 
    return pow(z.R*z.R + z.I*z.I, 0.5) ; 
} 

complex csquare(complex z)    //returns square of a complex 
{ 
    complex product ; 
    product.R = z.R*z.R - z.I*z.I ; 
    product.I = 2*z.R*z.I ; 
    return product ; 

} 

complex csum(complex z1, complex z2) // sums two complex numbers 
{ 
    complex sum ; 
    sum.R = z1.R + z2.R ; 
    sum.I = z1.I + z2.I; 
    return sum ; 
} 

Я получаю несколько реальных значений, они получают очень большой действительно быстро и лежат вне радиуса 2, затем много от реальных частей -nan и мнимых частей -nan. Любые предложения относительно того, чего я не хватает?

Answer 1

Я не уверен, что результаты вы ожидаете, но, возможно, это приборная выход поможет вам диагностировать проблему:

i = 0 
i = 0, j = 0 
i = 0, j = 0, k = 0: 0.000000 0.000000 
i = 0, j = 1 
i = 0, j = 1, k = 0: -2.100000 -1.200000 
i = 0, j = 2 
i = 0, j = 2, k = 0: 0.870000 3.853333 
i = 0, j = 3 
i = 0, j = 3, k = 0: -16.191278 5.531467 
i = 0, j = 4 
i = 0, j = 4, k = 0: 229.460353 -180.283027 
i = 0, j = 5 
i = 0, j = 5, k = 0: 20147.983719 -82736.760384 
i = 0, j = 6 
i = 0, j = 6, k = 0: -6439430272.999375 -3333957803.416253 
i = 0, j = 7 
i = 0, j = 7, k = 0: 30350987605860679680.000000 42937577616442236928.000000 
i = 0, j = 8 
i = 0, j = 8, k = 0: -922453122916892839668491877362832506880.000000 2606395772124638489891541615388582215680.000000 
i = 0, j = 9 
i = 0, j = 9, k = 0: -5942379156970062175257857153425511563624711669037998408592102022831643293646848.000000 -4808555839107517668369914051798230293111512606121703008400866410836391727464448.000000 
i = 0, j = 10 
i = 0, j = 10, k = 0: 12189660787377226979177299907109395027773453611835899247334853368367297050334425366457359535808690377953524893091770076537471791519164311649924318416907272192.000000 57148523986878398200502132726020983696472380197135570914789139569106551459309671849901109020634047615447493780748138450365838320365732961897986301575347306496.000000 
i = 0, j = 10, k = 1: nan inf 
i = 0, j = 10, k = 2: nan nan 

Обратите внимание, как внутренний цикл разбивается на первой итерации, все время , Я думаю, что проблема в том, что dz настроен на -2.1 - 1.2i для начала. Но без определения кривой Мандельброта, которую вы пытаетесь построить, я не могу сказать, что должно происходить.

Это было произведено в этом однофайльном варианте вашего кода, используя Complex вместо типа complex, который противоречит стандарту C99. Сложные функции являются статическими, потому что замолкает предупреждения компилятора при компиляции с:

$ gcc -O3 -g -std=c99 -Wall -Wextra -Wmissing-prototypes -Wstrict-prototypes cx.c -o cx 

Источник:

#include <stdio.h> 
#include <math.h> 

typedef struct Complex{ double R; double I; } Complex; 

static double cmag(Complex z) 
{ 
    return pow(z.R*z.R + z.I*z.I, 0.5); 
} 

static Complex csquare(Complex z)    //returns square of a complex 
{ 
    Complex product; 
    product.R = z.R*z.R - z.I*z.I; 
    product.I = 2*z.R*z.I; 
    return product; 
} 

static Complex csum(Complex z1, Complex z2) // sums two complex numbers 
{ 
    Complex sum; 
    sum.R = z1.R + z2.R; 
    sum.I = z1.I + z2.I; 
    return sum; 
} 

int main(void) 
{ 
    double dzrmin, dzrmax, dzimin, dzimax, dzr, dzi, r0, i0, maxrad; 
    int i, j, k, iter, Ndr, Ndi; 
    Complex z, dz; 
    dzrmin = -2.1; 
    dzrmax = 0.6; 
    dzimin = -1.2; 
    dzimax = 1.2; 
    Ndr = 200; 
    Ndi = 180; 
    dzr = (dzrmax - dzrmin)/Ndr; 
    dzi = (dzimax - dzimin)/Ndi; 
    r0 = 0.0; 
    i0 = 0.0; 
    z.R = r0; 
    z.I = i0; 
    dz.R = dzrmin; 
    dz.I = dzimin; 
    maxrad = 2.0; 
    iter = 80; 
    for (i = 0; i < Ndr; i++) 
    { 
     printf("i = %d\n", i); 
     for (j = 0; j < Ndi; j++) 
     { 
      printf("i = %d, j = %d\n", i, j); 
      for (k = 0; k < iter; k++) 
      { 
       printf("i = %d, j = %d, k = %d: ", i, j, k); 
       printf("%.6lf %.6lf\n", z.R, z.I); 
       z = csum(csquare(z), dz); 
       if (cmag(z) > maxrad) 
        break; 
      } 
      dz.I += dzi; 
     } 
     dz.R += dzr; 
     z.R = r0; 
     z.I = i0; 
    } 
    return 0; 
}