/* NAME * phase1d - plot phase-space of a one-dimensional map * NOTES * Includes the file "maps1d.c". * MAPS * The following four one-dimensional maps are allowed: * .IP Logistic\ Map: * x(t+1) = 4 * r * x(t) * (1.0 - x(t)) * .IP Tent\ Map: * x(t+1) = (x(t) <= 0.5) ? 2 * r * x(t) : 2r * (1.0 - x(t)) * .IP Sine\ Map: * x(t+1) = sin(x(t) * PI * aux * 2 * r) / 2 + 0.5 * .IP Gaussian\ Map: * x(t+1) = r * exp(-aux * (x(t) - 0.5) * (x(t) - 0.5)) * .LP * See the file "maps1d.c" to see how to add user-defined maps. * BUGS * No sanity checks are performed to make sure that any of the * options make sense. * AUTHOR * Copyright (c) 1997, Gary William Flake. * * Permission granted for any use according to the standard GNU * ``copyleft'' agreement provided that the author's comments are * neither modified nor removed. No warranty is given or implied. */ #include #include #include #include #include "misc.h" int width = 640, height = 480, points = 10, skip = 0; int mag = 1, invert = 0, arrows = 1; double r = 1.0, x0 = 0.123456, aux = 1.0, dx = 0.0; char *term = NULL, *func = "log"; char help_string[] = "\ A phase-space diagram is plotted for a one-dimensional map \ according to the specified options. If the option for -dx is \ non-zero, then two trajectories are plotted: one staring at x0, \ and the other starting at (x0 + dx). See the MAPS section of the \ manual page for details of what maps are supported. User defined \ maps can be added to the file maps1d.c, but you must recompile \ the program. \ "; OPTION options[] = { { "-width", OPT_INT, &width, "Width of the plot in pixels." }, { "-height", OPT_INT, &height, "Height of the plot in pixels." }, { "-points", OPT_INT, &points, "Number of points to plot." }, { "-skip", OPT_INT, &skip, "Number of points to skip." }, { "-r", OPT_DOUBLE, &r, "Value for r." }, { "-aux", OPT_DOUBLE, &aux, "Auxiliary map parameter." }, { "-x0", OPT_DOUBLE, &x0, "Initial value for x." }, { "-dx", OPT_DOUBLE, &dx, "Difference for second trajectory." }, { "-func", OPT_STRING, &func, "Map function to use (one of 'log', 'tent', 'sin', or 'gauss')." }, { "-arrows", OPT_SWITCH, &arrows, "Show arrows to indicate directions?" }, { "-inv", OPT_SWITCH, &invert, "Invert all colors?" }, { "-mag", OPT_INT, &mag, "Magnification factor." }, { "-term", OPT_STRING, &term, "How to plot points." }, { NULL, OPT_NULL, NULL, NULL } }; /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ /* This function plots a little triangular arrow head that is located midway between (x1, y1) and (x2, y2). */ void plot_arrow(double x1, double y1, double x2, double y2) { double xm, ym, off; /* 'off' is an offset for the points in the arrow head. This * should probably be a parameter or a constant. */ off = 0.01; /* Calculate the midpoint. */ xm = (x1 + x2) / 2.0; ym = (y1 + y2) / 2.0; /* Moving north. */ if(x1 == x2 && fabs(y1 - y2) > off && y2 > y1) { plot_line(xm, ym, xm - off, ym - off, 1); plot_line(xm, ym, xm + off, ym - off, 1); plot_line(xm - off, ym - off, xm + off, ym - off, 1); } /* Moving south. */ else if(x1 == x2 && fabs(y1 - y2) > off && y1 > y2) { plot_line(xm, ym, xm - off, ym + off, 1); plot_line(xm, ym, xm + off, ym + off, 1); plot_line(xm - off, ym + off, xm + off, ym + off, 1); } /* Moving east. */ else if(y1 == y2 && fabs(x1 - x2) > off && x2 > x1) { plot_line(xm, ym, xm - off, ym - off, 1); plot_line(xm, ym, xm - off, ym + off, 1); plot_line(xm - off, ym - off, xm - off, ym + off, 1); } /* Moving west. */ else if(y1 == y2 && fabs(x1 - x2) > off && x1 > x2) { plot_line(xm, ym, xm + off, ym - off, 1); plot_line(xm, ym, xm + off, ym + off, 1); plot_line(xm + off, ym - off, xm + off, ym + off, 1); } } /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ #include "maps1d.c" int main(int argc, char **argv) { extern int plot_mag; extern int plot_inverse; int i; double x, y, xt, yt, xinc, (*f)(double, double); get_options(argc, argv, options, help_string); plot_mag = mag; plot_inverse = invert; plot_init(width, height, 2, term); plot_set_range(0.0, 1.0, 0.0, 1.0); if(!term || strcmp(term, "ps") != 0) plot_set_all(0); f = get_named_function(func); /* First, draw the phase-space. */ xinc = 1.0 / (width - 1); xt = 0.0; yt = f(xt, r); for(i = 0, x = 0.0; i < width; i++, x += xinc) { y = f(x, r); plot_line(xt, yt, x, y, 1); xt = x; yt = y; } /* Draw the identity line. */ plot_line(0.0, 0.0, 1.0, 1.0, 1); /* Draw the trajectory. */ for(i = 0, x = x0; i < points + skip; i++, x = y) { y = f(x, r); if(i >= skip) { plot_line(x, x, x, y, 1); plot_line(x, y, y, y, 1); if(arrows) { plot_arrow(x, x, x, y); plot_arrow(x, y, y, y); } } } /* Draw the second trajectory, if appropriate. */ if(dx > 0.0) for(i = 0, x = x0 + dx; i < points + skip; i++, x = y) { y = f(x, r); if(i >= skip) { plot_line(x, x, x, y, 1); plot_line(x, y, y, y, 1); if(arrows) { plot_arrow(x, x, x, y); plot_arrow(x, y, y, y); } } } plot_finish(); exit(0); } /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */