Loading web-font TeX/Main/Regular

Visualizing \sqrt{2}, \pi, and e Using Random Walk 〜1 Million Digits〜 [C, gnuplot]

Saturday, October 20, 2018

C-language gnuplot Random Walk YouTube

t f B! P L

YouTube

 

Simulation [gnuplot]

Case 1 : \sqrt{2}

Make DAT file [C]

Get "r2_1M.dat" (line:21) → The Square Root of Two to 1 Million Digits (APOD)
▼ RW-r2.c
  1. #include <stdio.h>
  2. #include <math.h>
  3. #include <string.h>
  4. #include <stdlib.h>
  5. #define LIMIT 1000000 // Digits
  6. #define LENGTH 80 // The number of characters per line
  7. #define NUM0 4 // Number of leading zeros
  8. // Link(APOD) : 1.414213562373095048...
  9. // -> r2_1M.dat:0000414213562373095048...
  10. #define PI 3.1415
  11.  
  12. int main(void){
  13. FILE *fp, *fq;
  14. char str[256];
  15. double x = 0; // Position.x
  16. double y = 0; // Position.y
  17. int i, j, k; // For Loop counter
  18. int digit; // Current digit
  19. int num; // Number
  20. fp = fopen("r2_1M.dat", "r"); // Make this file by yourself
  21. fq = fopen("randomwalk-r2.dat", "w");
  22.  
  23. if(fp==NULL){
  24. printf("File not found. Program terminated.\n");
  25. exit(1);
  26. }
  27.  
  28. fprintf(fq, "%.3lf\t%.3lf\n",x, y); // Initiate value
  29. // i: row, j: column, k: 0 - 9
  30. for(i=1; fgets(str, sizeof(str), fp)!=NULL; i++){
  31. for (j=0; j<LENGTH; j++){
  32. // Calculate current digit
  33. digit = (LENGTH*(i-1)+j+1) - NUM0;
  34.  
  35. // Ignore leading zero NUM0 times
  36. if(digit<=0) continue;
  37. // From (NUM0+1)th digit on
  38. for(k=0; k<10; k++){ // Check matching between k and num
  39. num = str[j] - '0';
  40.  
  41. if(num ==k){ // Calculate x and y
  42. x += cos(2*PI*k/10);
  43. y += sin(2*PI*k/10);
  44. fprintf(fq, "%.3lf\t%.3lf\t%d\t%7d\n",x, y, k, digit-NUM0);
  45. break;
  46. }
  47. }
  48. if(digit==LIMIT) break;
  49. }
  50. if(digit==LIMIT) break;
  51. }
  52. fclose(fq);
  53. fclose(fp);
  54. return 0;
  55. }
→make "randomwalk-r2.dat"

Draw \sqrt{2} random walk [gnuplot]

▼ RW-r2.plt
  1. # Setting --------------------
  2. reset
  3. set terminal gif animate delay 6 size 1280, 720
  4. set output "randomwalk-pi.gif"
  5. set size ratio 1
  6. set nokey
  7.  
  8. set xr[-2000:200]
  9. set yr[-600:1600]
  10. set tics font 'Times, 20'
  11. unset xl
  12. unset yl
  13.  
  14. N = 1000000 # Digits
  15. dec = 200 # Decimation
  16.  
  17. #system "mkdir png" # Make folder for storing PNG files
  18.  
  19. DIGIT(i) = sprintf("Digits: %d\n", i) # Label displaying digits
  20.  
  21. # Plot --------------------
  22. do for[i=0:N]{
  23. # Digits
  24. set label 1 left DIGIT(i) font 'Times, 22' at graph 1.04, 0.03
  25. # Decimate and draw (the number of files: N/200+1)
  26. if(i%dec==0){
  27. # Draw initiate value
  28. if(i==0){
  29. plot 1/0
  30. continue
  31. }
  32. # Update
  33. plot 'randomwalk-r2.dat' every ::0::i-1 using 1:2 with lines lc rgb "red"
  34. }
  35. }
  36.  
  37. set out

▼ randomwalk-r2.gif (dec=5000)

Case 2 : \pi

Make DAT file [C]

Get "pi_1M.dat" (line:23) → 1 Million Digits of Pi (Pi Land by Eve)
▼ RW-pi.c
  1. #define LENGTH 50 // The number of characters per line
  2. #define NUM0 0 // Number of leading zeros
  3. // Link(Pi Land) :141592653589793238462643...
  4. // -> pi_1M.dat :141592653589793238462643...
  1. fp = fopen("pi_1M.dat", "r");
  2. fq = fopen("randomwalk-pi.dat", "w");
→make "randomwalk-pi.dat"

Draw \pi random walk [gnuplot]

▼ RW-pi.plt
  1. plot 'randomwalk-pi.dat' every ::0::i-1 using 1:2 with lines lc rgb "blue"

▼ randomwalk-pi.gif (dec=5000)

Case 3 : e

Make DAT file [C]

Get "e_1M.dat" (line:23) → The Number e to 1 Million Digits (APOD)
▼ RW-pi.c
  1. #define LENGTH 80 // The number of characters per line
  2. #define NUM0 4 // Number of leading zeros
  3. // Link(APOD) : 2.718281828459045235360...
  4. // -> e_1M.dat :0000718281828459045235360...
  1. fp = fopen("e_1M.dat", "r");
  2. fq = fopen("randomwalk-e.dat", "w");
→make "randomwalk-e.dat"

Draw e random walk [gnuplot]

▼ RW-e.plt
  1. plot 'randomwalk-e.dat' every ::0::i-1 using 1:2 with lines lc rgb "green"

▼ randomwalk-e.gif (dec=5000)

Extra

Draw random walks of \sqrt{2}, \pi, and e simultaneously

Read DAT files created already and plot them

▼ 3RWs.plt
  1. # Setting --------------------
  2. reset
  3. set terminal png enhanced font "Times" 20 size 720, 720
  4. set xr[-2000:200]
  5. set yr[-800:1400]
  6. set size ratio 1
  7. set tics font 'Times, 20'
  8. set xtics 400
  9. set ytics 400
  10. unset key
  11.  
  12. # Axes
  13. set arrow 1 nohead from 0, -800 to 0, 1400 lc rgb "black" lw 2
  14. set arrow 2 nohead from -2000, 0 to 200, 0 lc rgb "black" lw 2
  15.  
  16. # Digits
  17. N = 1000000
  18.  
  19. # Arrays
  20. array num[3] = ["r2", "pi", "e"] # Number
  21. array color[3] = ["red", "blue", "green"] # Color
  22.  
  23. # Functions --------------------
  24. PNG(i) = sprintf("img_%04d.png", i) # Output img_000i.png
  25. DAT(i) = sprintf("randomwalk-%s.dat", num[i]) # num[i] -> name of DAT files
  26.  
  27. # Plot DAT files --------------------
  28. do for[i=1:3]{
  29. set output PNG(i)
  30. str = "plot " # Command (need space)
  31.  
  32. # Make the command using string concatenation ("A"."B" -> "AB")
  33. do for [j=1:i] {
  34. if(j!=1){
  35. str = str.", " # Need ", " for plotting DAT files together
  36. }
  37. str = str.sprintf("'%s' every ::0::N with lines lc rgb '%s' lw 0.5", DAT(j), color[j])
  38. }
  39.  
  40. eval str # Execute this command
  41. # i=1: plot 'randomwalk-r2.dat' every ::0::N with lines lc rgb 'red' lw 0.5
  42. # i=2: plot 'randomwalk-r2.dat' every ::0::N with lines lc rgb 'red' lw 0.5,
  43. # 'randomwalk-pi.dat' every ::0::N with lines lc rgb 'blue' lw 0.5
  44. # i=3: plot 'randomwalk-r2.dat' every ::0::N with lines lc rgb 'red' lw 0.5,
  45. # 'randomwalk-pi.dat' every ::0::N with lines lc rgb 'blue' lw 0.5,
  46. # 'randomwalk-e.dat' every ::0::N with lines lc rgb 'green' lw 0.5
  47. }
  48.  
  49. set out

Output PNG files

img_0001.png (\sqrt{2}) img_0002.png (\sqrt{2}, \pi) img_0003.png (\sqrt{2}, \pi, e)

Search This Blog

Translate

QooQ