Advertisement
Magog

PHP bug

Nov 29th, 2012
116
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
PHP 17.64 KB | None | 0 0
  1. <?php
  2. $text = "== {{int:filedesc}} ==\r\n{{Information\r\n|Description    ={{en|1=Parabolic Julia set for internal angle 1 over 5 }}\r\n|Source         =own work with help of : see comments in code \r\n|Author         =[[User:Adam majewski|Adam majewski]]\r\n|Date           =2012-10-07\r\n|Permission     =\r\n|other_versions =\r\n}}\r\n\r\n==C src code ==\r\n<source lang=c>\r\n/*\r\n\r\nAdam Majewski\r\nfraktal.republika.pl\r\n\r\nc console progam using \r\n* symmetry\r\n* openMP\r\n\r\nIt uses modified DEM method to draw parabolic julia sets\r\n\r\n\r\n\r\ngcc t.c -lm -Wall -fopenmp -march=native \r\ntime ./a.out\r\n\r\n\r\nFile s1000000f1.5.pgm saved. \r\nCx  = 0.356763 \r\nCy  = 0.328582 \r\nalfax  = 0.154508 \r\nalfay  = 0.475528 \r\ndistorsion of image  = 1.000000 \r\n\r\nreal\t123m19.106s\r\n----------------------------\r\nFile s10000000f1.5.pgm saved. \r\nCx  = 0.356763 \r\nCy  = 0.328582 \r\nalfax  = 0.154508 \r\nalfay  = 0.475528 \r\ndistorsion of image  = 1.000000 \r\n\r\nreal\t1023m48.596s = 17 fours\r\n-------------------------------------\r\nso s100000000f1.5.pgm shoud take 17 days !!!!!! \r\n\r\n\r\n\r\n\r\n\r\n\r\n*/\r\n\r\n\r\n\r\n#include <stdio.h>\r\n#include <stdlib.h> // malloc\r\n#include <string.h> // strcat\r\n#include <math.h> // M_PI; needs -lm also \r\n#include <complex.h>\r\n#include <omp.h> // OpenMP; needs also -fopenmp\r\n\r\n\r\n/* --------------------------------- global variables and consts ------------------------------------------------------------ */\r\n\r\n\r\n// virtual 2D array and integer ( screen) coordinate\r\n// Indexes of array starts from 0 not 1 \r\nunsigned int ix, iy; // var\r\nunsigned int ixMin = 0; // Indexes of array starts from 0 not 1\r\nunsigned int ixMax ; //\r\nunsigned int iWidth ; // horizontal dimension of array\r\nunsigned int ixAxisOfSymmetry  ; // \r\nunsigned int iyMin = 0; // Indexes of array starts from 0 not 1\r\nunsigned int iyMax ; //\r\nunsigned int iyAxisOfSymmetry  ; // \r\nunsigned int iyAbove ; // var, measured from 1 to (iyAboveAxisLength -1)\r\nunsigned int iyAboveMin = 1 ; //\r\nunsigned int iyAboveMax ; //\r\nunsigned int iyAboveAxisLength ; //\r\nunsigned int iyBelowAxisLength ; //\r\nunsigned int iHeight = 1000; //  odd number !!!!!! = (iyMax -iyMin + 1) = iyAboveAxisLength + iyBelowAxisLength +1\r\n// The size of array has to be a positive constant integer \r\nunsigned int iSize ; // = iWidth*iHeight; \r\n\r\n\r\n// memmory 1D array \r\nunsigned char *data;\r\n// unsigned int i; // var = index of 1D array\r\nunsigned int iMin = 0; // Indexes of array starts from 0 not 1\r\nunsigned int iMax ; // = i2Dsize-1  = \r\n// The size of array has to be a positive constant integer \r\n// unsigned int i1Dsize ; // = i2Dsize  = (iMax -iMin + 1) =  ;  1D array with the same size as 2D array\r\n\r\n\r\n/* world ( double) coordinate = dynamic plane */\r\n  const double ZxMin=-1.5;\r\n  const double ZxMax=1.5;\r\n  const double ZyMin=-1.5;\r\n  const double ZyMax=1.5;\r\n  double PixelWidth; // =(ZxMax-ZxMin)/iXmax;\r\n  double PixelHeight; // =(ZyMax-ZyMin)/iYmax;\r\n  double distanceMax;\r\n  double ratio ;\r\n  double lambda=1.5;\r\n\r\n\r\n\r\n// complex numbers of parametr plane \r\ndouble Cx; // c =Cx +Cy * i\r\ndouble Cy;\r\ndouble complex c; // \r\n\r\ndouble complex alfa; // alfa fixed point alfa=f(alfa)\r\n\r\nunsigned long int iterMax  = 10000000; //iHeight*100;\r\ndouble ER = 2.0; // Escape Radius for bailout test \r\ndouble ER2;\r\n\r\n\r\n/* ------------------------------------------ functions -------------------------------------------------------------*/\r\n\r\n\r\n\r\n/* find c in component of Mandelbrot set \r\n \r\n uses code by Wolf Jung from program Mandel\r\n see function mndlbrot::bifurcate from mandelbrot.cpp\r\n http://www.mndynamics.com/indexp.html\r\n\r\n  */\r\ndouble complex GiveC(double InternalAngleInTurns, double InternalRadius, unsigned int period)\r\n{\r\n  //0 <= InternalRay<= 1\r\n  //0 <= InternalAngleInTurns <=1\r\n  double t = InternalAngleInTurns *2*M_PI; // from turns to radians\r\n  double R2 = InternalRadius * InternalRadius;\r\n  //double Cx, Cy; /* C = Cx+Cy*i */\r\n  switch ( period ) // of component \r\n  {\r\n    case 1: // main cardioid\r\n      Cx = (cos(t)*InternalRadius)/2-(cos(2*t)*R2)/4; \r\n      Cy = (sin(t)*InternalRadius)/2-(sin(2*t)*R2)/4; \r\n      break;\r\n   case 2: // only one component \r\n      Cx = InternalRadius * 0.25*cos(t) - 1.0;\r\n      Cy = InternalRadius * 0.25*sin(t); \r\n      break;\r\n  // for each period  there are 2^(period-1) roots. \r\n  default: // higher periods : to do\r\n      Cx = 0.0;\r\n      Cy = 0.0; \r\n    break; }\r\n\r\n  return Cx + Cy*I;\r\n}\r\n\r\n\r\n/*\r\n\r\nhttp://en.wikipedia.org/wiki/Periodic_points_of_complex_quadratic_mappings\r\nz^2 + c = z\r\nz^2 - z + c = 0\r\nax^2 +bx + c =0 // ge3neral for  of quadratic equation\r\nso :\r\na=1\r\nb =-1\r\nc = c\r\nso :\r\n\r\nThe discriminant is the  d=b^2- 4ac \r\n\r\nd=1-4c = dx+dy*i\r\nr(d)=sqrt(dx^2 + dy^2)\r\nsqrt(d) = sqrt((r+dx)/2)+-sqrt((r-dx)/2)*i = sx +- sy*i\r\n\r\nx1=(1+sqrt(d))/2 = beta = (1+sx+sy*i)/2\r\n\r\nx2=(1-sqrt(d))/2 = alfa = (1-sx -sy*i)/2\r\n\r\nalfa : attracting when c is in main cardioid of Mandelbrot set, then it is in interior of Filled-in Julia set, \r\nit means belongs to Fatou set ( strictly to basin of attraction of finite fixed point )\r\n\r\n*/\r\n// uses global variables : \r\n//  ax, ay (output = alfa(c)) \r\ndouble complex GiveAlfaFixedPoint(double complex c)\r\n{\r\n  double dx, dy; //The discriminant is the  d=b^2- 4ac = dx+dy*i\r\n  double r; // r(d)=sqrt(dx^2 + dy^2)\r\n  double sx, sy; // s = sqrt(d) = sqrt((r+dx)/2)+-sqrt((r-dx)/2)*i = sx + sy*i\r\n  double ax, ay;\r\n \r\n // d=1-4c = dx+dy*i\r\n dx = 1 - 4*creal(c);\r\n dy = -4 * cimag(c);\r\n // r(d)=sqrt(dx^2 + dy^2)\r\n r = sqrt(dx*dx + dy*dy);\r\n //sqrt(d) = s =sx +sy*i\r\n sx = sqrt((r+dx)/2);\r\n sy = sqrt((r-dx)/2);\r\n // alfa = ax +ay*i = (1-sqrt(d))/2 = (1-sx + sy*i)/2\r\n ax = 0.5 - sx/2.0;\r\n ay =  sy/2.0;\r\n \r\n\r\nreturn ax+ay*I;\r\n}\r\n\r\n\r\ndouble GiveDistanceBetween(double complex z1, double complex z2 )\r\n{double dx,dy;\r\n \r\n dx = creal(z1) - creal(z2);\r\n dy = cimag(z1) - cimag(z2);\r\n return sqrt(dx*dx+dy*dy);\r\n \r\n} \r\n\r\n\r\n\r\nint setup()\r\n{\r\n\r\n  unsigned int period = 5; // period of secondary component joined by root point\r\n  unsigned int denominator;\r\n  double InternalAngle;\r\n  \r\n\r\n  denominator = period;\r\n  InternalAngle = 1.0/((double) denominator);\r\n\r\n  c = GiveC(InternalAngle, 1.0, 1) ;\r\n  Cx=creal(c);\r\n  Cy=cimag(c);\r\n  alfa = GiveAlfaFixedPoint(c);\r\n\r\n  /* 2D array ranges */\r\n  if (!(iHeight % 2)) iHeight+=1; // it sholud be even number (variable % 2) or (variable & 1)\r\n  iWidth = iHeight;\r\n  iSize = iWidth*iHeight; // size = number of points in array \r\n  // iy\r\n  iyMax = iHeight - 1 ; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].\r\n  iyAboveAxisLength = (iHeight -1)/2;\r\n  iyAboveMax = iyAboveAxisLength ; \r\n  iyBelowAxisLength = iyAboveAxisLength; // the same \r\n  iyAxisOfSymmetry = iyMin + iyBelowAxisLength ; \r\n  // ix\r\n  \r\n  ixMax = iWidth - 1;\r\n\r\n  /* 1D array ranges */\r\n  // i1Dsize = i2Dsize; // 1D array with the same size as 2D array\r\n  iMax = iSize-1; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].\r\n\r\n\r\n /* Pixel sizes */\r\n  PixelWidth = (ZxMax-ZxMin)/ixMax; //  ixMax = (iWidth-1)  step between pixels in world coordinate \r\n  PixelHeight = (ZyMax-ZyMin)/iyMax;\r\n  ratio = ((ZxMax-ZxMin)/(ZyMax-ZyMin))/((float)iWidth/(float)iHeight); // it should be 1.000 ...\r\n  distanceMax = PixelWidth; \r\n  //\r\n  \r\n  ER2 = ER * ER;\r\n\r\n  /* create dynamic 1D arrays for colors ( shades of gray ) */\r\n  \r\n  data = malloc( iSize * sizeof(unsigned char) );\r\n  if (data == NULL )\r\n    {\r\n      fprintf(stderr,\" Could not allocate memory\\n\");\r\n      getchar();\r\n      return 1;\r\n    }\r\n  else fprintf(stderr,\" memory is OK \\n\");\r\n\r\n  \r\n \r\n  return 0;\r\n\r\n}\r\n\r\n\r\n\r\n// from screen to world coordinate ; linear mapping\r\n// uses global cons\r\ndouble GiveZx(unsigned int ix)\r\n{ return (ZxMin + ix*PixelWidth );}\r\n\r\n// uses globaal cons\r\ndouble GiveZy(unsigned int iy)\r\n{ return (ZyMax - iy*PixelHeight);} // reverse y axis\r\n\r\n\r\n// forward iteration of complex quadratic polynomial\r\n// fc(z) = z*z +c\r\n// z0 = critical point \r\n// uses global var \r\nlong int GiveLastIteration(double Zx, double Zy )\r\n{ \r\n    long int iter;\r\n    // z= Zx + Zy*i\r\n    double Zx2, Zy2; /* Zx2=Zx*Zx;  Zy2=Zy*Zy  */\r\n \r\n     /* initial value of orbit = critical point Z= 0 */\r\n                       \r\n                        Zx2=Zx*Zx;\r\n                        Zy2=Zy*Zy;\r\n                        // for iter from 0 to iterMax because of ++ after last loop\r\n                        for (iter=0; iter<iterMax && ((Zx2+Zy2)<ER2); ++iter)\r\n                        {\r\n                            Zy=2*Zx*Zy + Cy;\r\n                            Zx=Zx2-Zy2 + Cx;\r\n                            Zx2=Zx*Zx;\r\n                            Zy2=Zy*Zy;\r\n                        };\r\n  \r\n   return iter;\r\n}\r\n\r\n\r\n\r\n\r\n \r\n\r\n\r\n\r\n/*\r\n estimates distance from point c to nearest point in Julia  set \r\n for Fc(z)= z*z + c\r\n z(n+1) = Fc(zn)  \r\n this function is based on function  mndlbrot::dist  from  mndlbrot.cpp\r\n from program mandel by Wolf Jung (GNU GPL )\r\n http://www.mndynamics.com/indexp.html \r\n\r\nHyunsuk Kim  : \r\nFor Julia sets, z is the variable and c is a constant. Therefore df[n+1](z)/dz = 2*f[n]*f'[n] -- you don't add 1.\r\n-----------------\r\n\r\n\r\n\"The following method was taught to me by Christian Henriksen.\r\nIt is a distance estimator method (different from Milnor's).\r\nWhen iterating forward, keep also tract of the value of the derivative.\r\nMore precisely when computing z_n from z_0, also compute dz_n/dz_0.\r\nby the chain rule, it is equal to the product of the derivatives of f taken at z_0, z_1, .. , up to z_n-1\r\nTherefore you can compute it recursively alongside and like the value of z_n.\r\nThen choose a point in the Julia set, here I chose the critical point.\r\nThen if the ball :\r\n- of center z_n \r\n- and radius lambda * size of a pixel * derivative \r\ncontains this point there is a good chance \r\nthat the pixel of center z_0 \r\ncontains an iterated preimage of the critical point, thus I color z_0 in black.\r\n\r\nHere lambda is a thickness factor that you must choose. \r\nnot too small otherwise you do not see enough points, but not too big because otherwise you get artifacts.\r\n\r\nBest Regards,\r\nArnaud.\"\r\n\r\n\r\n */\r\n int jdist(double Zx, double Zy, double Cx, double Cy ,  int iter_max)\r\n { \r\n int i;\r\n double x = Zx, /* Z = x+y*i */\r\n         y = Zy, \r\n         /* Zp = xp+yp*1 = 1  */\r\n         xp = 1, \r\n         yp = 0, \r\n         /* temporary */\r\n         nz,  \r\n         nzp,\r\n         /* a = abs(z) */\r\n         //a,\r\n         radius,\r\n         distance; \r\n\r\n\r\n \r\n for (i = 1; i <= iter_max; i++)\r\n  { /* first derivative   zp = 2*z*zp  = xp + yp*i; */\r\n    \r\n    nz = 2*(x*xp - y*yp) ; \r\n    yp = 2*(x*yp + y*xp); \r\n    xp = nz;\r\n    /* z = z*z + c = x+y*i */\r\n    nz = x*x - y*y + Cx; \r\n    y = 2*x*y + Cy; \r\n    x = nz; \r\n    //\r\n    nz = (x*x + y*y); // abs2(z)\r\n    nzp = (xp*xp + yp*yp);\r\n    if ( nz > ER2) return 0; // if escapes //nzp > 1e60 ||\r\n    if (nzp>1e60) return 0; // ? interior ?\r\n    distance=GiveDistanceBetween(alfa,x+y*I);\r\n    nzp = sqrt(xp*xp + yp*yp);\r\n    radius=lambda*PixelWidth*nzp;\r\n    \r\n    if (distance<radius)return 1; // \r\n  }\r\n \r\n return 0; // \r\n }\r\n\r\n\r\nunsigned char GiveColor(unsigned int ix, unsigned int iy)\r\n{ \r\n   double Zx, Zy; //  Z= Zx+ZY*i;\r\n  unsigned char color; // gray from 0 to 255 \r\n// unsigned char ColorList[]={255,230,180}; /* shades of gray used in imaage\r\n  //color=0;\r\n  // from screen to world coordinate \r\n  Zx = GiveZx(ix);\r\n  Zy = GiveZy(iy);\r\n  //if ( GiveLastIteration(Zx, Zy ) == iterMax)\r\n    // color = 0; // interior\r\n    // else { \r\n\r\n// only dist\r\n  if (jdist(Zx, Zy, Cx, Cy , iterMax))\r\n                 color = 0; // black = boundary = Julia set\r\n                 else color = 255; // Fatou set\r\n          //}\r\n  return color;\r\n}\r\n\r\n\r\n/* -----------  array functions -------------- */\r\n\r\n\r\n/* gives position of 2D point (iX,iY) in 1D array  ; uses also global variable iWidth */\r\nunsigned int Give_i(unsigned int ix, unsigned int iy)\r\n{ return ix + iy*iWidth; }\r\n//  ix = i % iWidth;\r\n//  iy = (i- ix) / iWidth;\r\n//  i  = Give_i(ix, iy);\r\n\r\n\r\n\r\n\r\n// plots raster point (ix,iy) \r\nint PlotPoint(unsigned int ix, unsigned int iy, unsigned char iColor)\r\n{\r\n unsigned i; /* index of 1D array */\r\n i = Give_i(ix,iy); /* compute index of 1D array from indices of 2D array */\r\n data[i] = iColor;\r\n\r\nreturn 0;\r\n}\r\n\r\n\r\n// fill array \r\n// uses global var :  ...\r\n// scanning complex plane \r\nint FillArray(unsigned char data[] )\r\n{\r\n  unsigned int ix, iy; // pixel coordinate \r\n\r\n\r\n// for all pixels of image \r\nfor(iy = iyMin; iy<=iyMax; ++iy) \r\n  { printf(\" %d z %d\\n\", iy, iyMax); //info \r\n    for(ix= ixMin; ix<=ixMax; ++ix) PlotPoint(ix, iy, GiveColor(ix, iy) ); //  \r\n   } \r\n   \r\n return 0;\r\n}\r\n\r\n\r\n// fill array using symmetry of image \r\n// uses global var :  ...\r\nint FillArraySymmetric(unsigned char data[] )\r\n{\r\n   \r\n unsigned char Color; // gray from 0 to 255 \r\n\r\nprintf(\"axis of symmetry \\n\"); \r\niy = iyAxisOfSymmetry; \r\n#pragma omp parallel for schedule(dynamic) private(ix,Color) shared(ixMin,ixMax, iyAxisOfSymmetry)\r\nfor(ix=ixMin;ix<=ixMax;++ix) {//printf(\" %d from %d\\n\", ix, ixMax); //info  \r\n                              PlotPoint(ix, iy, GiveColor(ix, iy));\r\n}\r\n\r\n\r\n/*\r\nThe use of ‘shared(variable, variable2) specifies that these variables should be shared among all the threads.\r\nThe use of ‘private(variable, variable2)’ specifies that these variables should have a seperate instance in each thread.\r\n*/\r\n\r\n#pragma omp parallel for schedule(dynamic) private(iyAbove,ix,iy,Color) shared(iyAboveMin, iyAboveMax,ixMin,ixMax, iyAxisOfSymmetry)\r\n\r\n// above and below axis \r\nfor(iyAbove = iyAboveMin; iyAbove<=iyAboveMax; ++iyAbove) \r\n  {printf(\" %d from %d\\n\", iyAbove, iyAboveMax); //info \r\n  for(ix=ixMin; ix<=ixMax; ++ix) \r\n\r\n  { // above axis compute color and save it to the array\r\n    iy = iyAxisOfSymmetry + iyAbove;\r\n    Color = GiveColor(ix, iy);\r\n    PlotPoint(ix, iy, Color ); \r\n    // below the axis only copy Color the same as above without computing it \r\n    PlotPoint(ixMax-ix, iyAxisOfSymmetry - iyAbove , Color ); \r\n   } \r\n}  \r\n return 0;\r\n}\r\n\r\n\r\n\r\n\r\n// Check Orientation of image : first quadrant in upper right position\r\n// uses global var :  ...\r\nint CheckOrientation(unsigned char data[] )\r\n{\r\n   unsigned int ix, iy; // pixel coordinate \r\n   double Zx, Zy; //  Z= Zx+ZY*i;\r\n   unsigned i; /* index of 1D array */\r\n   for(iy=iyMin;iy<=iyMax;++iy) \r\n   {\r\n     Zy = GiveZy(iy);\r\n    for(ix=ixMin;ix<=ixMax;++ix) \r\n   {\r\n\r\n    // from screen to world coordinate \r\n    Zx = GiveZx(ix);\r\n     i = Give_i(ix, iy); /* compute index of 1D array from indices of 2D array */\r\n     if (Zx>0 && Zy>0) data[i]=255-data[i];   // check the orientation of Z-plane by marking first quadrant */\r\n\r\n    }\r\n   }\r\n   \r\n  return 0;\r\n}\r\n\r\n\r\n\r\n// save data array to pgm file \r\nint SaveArray2PGMFile( unsigned char data[])\r\n{\r\n  \r\n  FILE * fp;\r\n  const unsigned int MaxColorComponentValue=255; /* color component is coded from 0 to 255 ;  it is 8 bit color file */\r\n  char name [10]; /* name of file */\r\n  sprintf(name,\"s%luf%3.1f\", iterMax, lambda); /*  */\r\n  char *filename =strcat(name,\".pgm\");\r\n  char *comment=\"# \";/* comment should start with # */\r\n\r\n  /* save image to the pgm file  */      \r\n  fp= fopen(filename,\"wb\"); /*create new file,give it a name and open it in binary mode  */\r\n  fprintf(fp,\"P5\\n %s\\n %u %u\\n %u\\n\",comment,iWidth,iHeight,MaxColorComponentValue);  /*write header to the file*/\r\n  fwrite(data,iSize,1,fp);  /*write image data bytes to the file in one step */\r\n  printf(\"File %s saved. \\n\", filename);\r\n  fclose(fp);\r\n\r\n  return 0;\r\n}\r\n\r\n\r\nint info()\r\n{\r\n // diplay info messages\r\n  printf(\"Cx  = %f \\n\", Cx); \r\n  printf(\"Cy  = %f \\n\", Cy);\r\n  // \r\n \r\n  printf(\"alfax  = %f \\n\", creal(alfa));\r\n  printf(\"alfay  = %f \\n\", cimag(alfa));\r\n  printf(\"distorsion of image  = %f \\n\", ratio);\r\nreturn 0;\r\n}\r\n\r\n\r\n/* -----------------------------------------  main   -------------------------------------------------------------*/\r\nint main()\r\n{\r\n// here are procedures for creating image file\r\n\r\n  setup();\r\n  // compute colors of pixels = image\r\n  //FillArray( data ); // no symmetry\r\n  FillArraySymmetric(data); \r\n  //  CheckOrientation( data );\r\n  SaveArray2PGMFile( data); // save array (image) to pgm file \r\n  free(data);\r\n  info();\r\n  return 0;\r\n}\r\n</source>\r\n\r\n[[Category:Uploaded with UploadWizard]]\r\n[[Category:Parabolic Julia sets]]\r\n[[Category:Images with C source code]]\r\n[[Category:Black and white]]\r\n[[Category:Complex quadratic map]]\r\n\r\n== {{int:license-header}} ==\r\n{{self|cc-by-sa-3.0|GFDL}}.";
  3.  
  4. echo "\r\n\r\n***Outputting the text in order to show that this is a valid string***\r\n";
  5. echo $text;
  6.  
  7. echo "\r\n\r\n***Outputting a regular expression run on the text, which returns null***\r\n";
  8. var_dump(preg_replace("/([\s\S]+?\s*?)\s*\s*(\{\{)/u", "Hello world", $text));
  9. ?>
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement