Line data Source code
1 : /* Copyright (C) 2016 The PARI group.
2 :
3 : This file is part of the PARI/GP package.
4 :
5 : PARI/GP is free software; you can redistribute it and/or modify it under the
6 : terms of the GNU General Public License as published by the Free Software
7 : Foundation; either version 2 of the License, or (at your option) any later
8 : version. It is distributed in the hope that it will be useful, but WITHOUT
9 : ANY WARRANTY WHATSOEVER.
10 :
11 : Check the License for details. You should have received a copy of it, along
12 : with the package; see the file 'COPYING'. If not, write to the Free Software
13 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14 :
15 : #include "pari.h"
16 : #include "paripriv.h"
17 :
18 : /********************************************************************/
19 : /** **/
20 : /** LOW-RES PLOT **/
21 : /** **/
22 : /********************************************************************/
23 : #define ISCR 64
24 : #define JSCR 22
25 :
26 : INLINE long
27 704 : DTOL(double t) { return (long)(t + 0.5); }
28 :
29 : static char
30 704 : PICT(long j) {
31 704 : switch(j%3) {
32 211 : case 0: return '_';
33 222 : case 1: return 'x';
34 271 : default: return '"';
35 : }
36 : }
37 : static char
38 11 : PICTZERO(long j) {
39 11 : switch(j%3) {
40 0 : case 0: return ',';
41 5 : case 1: return '-';
42 6 : default: return '`';
43 : }
44 : }
45 :
46 : static char *
47 22 : dsprintf9(double d, char *buf)
48 : {
49 22 : int i = 10;
50 :
51 47 : while (--i >= 0) {
52 47 : sprintf(buf, "%9.*g", i, d);
53 47 : if (strlen(buf) <= 9) break;
54 : }
55 22 : return buf;
56 : }
57 :
58 : typedef unsigned char screen[ISCR+1][JSCR+1];
59 :
60 : static void
61 693 : fill_gap(screen scr, long i, int jnew, int jpre)
62 : {
63 : int mid, i_up, i_lo, up, lo;
64 :
65 693 : if (jpre < jnew - 2) {
66 40 : up = jnew - 1; i_up = i;
67 40 : lo = jpre + 1; i_lo = i - 1;
68 653 : } else if (jnew < jpre - 2) {
69 60 : up = jpre - 1; i_up = i - 1;
70 60 : lo = jnew + 1; i_lo = i;
71 593 : } else return; /* if gap < 2, leave it as it is. */
72 :
73 100 : mid = (jpre+jnew)/2;
74 100 : if (mid>JSCR) mid=JSCR; else if (mid<0) mid=0;
75 100 : if (lo<0) lo=0;
76 250 : if (lo<=JSCR) while (lo <= mid) scr[i_lo][lo++] = ':';
77 100 : if (up>JSCR) up=JSCR;
78 240 : if (up>=0) while (up > mid) scr[i_up][up--] = ':';
79 : }
80 :
81 : static double
82 22 : todbl(GEN x) { return rtodbl(gtofp(x, LOWDEFAULTPREC)); }
83 :
84 : void
85 11 : pariplot(void* E, GEN (*fun)(void *E, GEN x), GEN a, GEN b, GEN ysmlu,GEN ybigu, long prec)
86 : {
87 11 : const char BLANK = ' ', YY = '|', XX_UPPER = '\'', XX_LOWER = '.';
88 : long jz, j, i, sig;
89 11 : pari_sp av = avma;
90 11 : int jnew, jpre = 0; /* for lint */
91 : GEN x, dx;
92 : double diff, dyj, ysml, ybig, y[ISCR+1];
93 : screen scr;
94 : char buf[80], z;
95 :
96 11 : sig=gcmp(b,a); if (!sig) return;
97 11 : if (sig<0) { x=a; a=b; b=x; }
98 11 : x = gtofp(a, prec);
99 11 : dx = divru(gtofp(gsub(b,a),prec), ISCR-1);
100 253 : for (j=1; j<=JSCR; j++) scr[1][j]=scr[ISCR][j]=YY;
101 693 : for (i=2; i<ISCR; i++)
102 : {
103 682 : scr[i][1] = XX_LOWER;
104 682 : scr[i][JSCR]= XX_UPPER;
105 14322 : for (j=2; j<JSCR; j++) scr[i][j] = BLANK;
106 : }
107 11 : ysml = ybig = 0.; /* -Wall */
108 715 : for (i=1; i<=ISCR; i++)
109 : {
110 704 : pari_sp av2 = avma;
111 704 : y[i] = gtodouble( fun(E, x) );
112 704 : set_avma(av2);
113 704 : if (i == 1)
114 11 : ysml = ybig = y[1];
115 : else
116 : {
117 693 : if (y[i] < ysml) ysml = y[i];
118 693 : if (y[i] > ybig) ybig = y[i];
119 : }
120 704 : x = addrr(x,dx);
121 : }
122 11 : set_avma(av);
123 11 : if (ysmlu) ysml = gtodouble(ysmlu);
124 11 : if (ybigu) ybig = gtodouble(ybigu);
125 11 : diff = ybig - ysml;
126 11 : if (!diff) { ybig += 1; diff= 1.; }
127 11 : dyj = ((JSCR-1)*3+2) / diff;
128 : /* work around bug in gcc-4.8 (32bit): plot(x=-5,5,sin(x)))) */
129 11 : jz = 3 - (long)(ysml*dyj + 0.5); /* 3 - DTOL(ysml*dyj) */
130 11 : z = PICTZERO(jz); jz /= 3;
131 715 : for (i=1; i<=ISCR; i++)
132 : {
133 704 : if (0<=jz && jz<=JSCR) scr[i][jz]=z;
134 704 : j = 3 + DTOL((y[i]-ysml)*dyj);
135 704 : jnew = j/3;
136 704 : if (i > 1) fill_gap(scr, i, jnew, jpre);
137 704 : if (0<=jnew && jnew<=JSCR) scr[i][jnew] = PICT(j);
138 704 : jpre = jnew;
139 : }
140 11 : pari_putc('\n');
141 11 : pari_printf("%s ", dsprintf9(ybig, buf));
142 715 : for (i=1; i<=ISCR; i++) pari_putc(scr[i][JSCR]);
143 11 : pari_putc('\n');
144 231 : for (j=(JSCR-1); j>=2; j--)
145 : {
146 220 : pari_puts(" ");
147 14300 : for (i=1; i<=ISCR; i++) pari_putc(scr[i][j]);
148 220 : pari_putc('\n');
149 : }
150 11 : pari_printf("%s ", dsprintf9(ysml, buf));
151 715 : for (i=1; i<=ISCR; i++) pari_putc(scr[i][1]);
152 11 : pari_putc('\n');
153 : {
154 : char line[10 + 32 + 32 + ISCR - 9];
155 11 : sprintf(line, "%10s%-9.7g%*.7g\n"," ",todbl(a),ISCR-9,todbl(b));
156 11 : pari_printf(line);
157 : }
158 : }
159 :
160 : void
161 11 : pariplot0(GEN a, GEN b, GEN code, GEN ysmlu,GEN ybigu, long prec)
162 : {
163 11 : push_lex(gen_0, code);
164 11 : pariplot((void*)code, &gp_eval, a, b, ysmlu, ybigu, prec);
165 11 : pop_lex(1);
166 11 : }
|
