#include <chplot.h> #include <math.h> #define NVAR 4 #define POINTS 50 int main() { void derivs(double t, double y[], double dydt[]) { dydt[0] = -y[1]; dydt[1]=y[0]-(1.0/t)*y[1]; dydt[2]=y[1]-(2.0/t)*y[2]; dydt[3]=y[2]-(3.0/t)*y[3]; } double t0=1, tf=10, y0[NVAR]; double t[POINTS], y1[NVAR][POINTS]; int points; int datasetnum=0, line_type=4, line_width = 2; array double theta2[360], r2[360]; array double x3[360], y3[360], z3[2][360]; double theta4[36], z4[20], r4[720]; double r, x5[30], y5[30], z5[900]; double x6[40], y6[60], z6[2400]; int i, j; class CPlot subplot, *spl; /* plot 1 */ y0[0]=j0(t0); y0[1]=j1(t0); y0[2]=jn(2,t0); y0[3]=jn(3,t0); points = odesolve(t, y1, derivs, t0, tf, y0); /* plot 2 */ linspace(theta2, 0, M_PI); r2 = sin(5*theta2); /* plot 3 */ linspace(x3, 0, 360); y3 = sin(x3*M_PI/180); for(i=0; i<360; i++) { z3[0][i] = cos(x3[i]*M_PI/180); z3[1][i] = y3[i]; } /* plot 4 */ linspace(theta4, 0, 360); linspace(z4, 0, 2*M_PI); for(i=0; i<36; i++) { for(j=0; j<20; j++) { r4[i*20+j] = 2+cos(z4[j]); } } /* plot 5 */ linspace(x5, -10, 10); linspace(y5, -10, 10); for(i=0; i<30; i++) { for(j=0; j<30; j++) { r = sqrt(x5[i]*x5[i]+y5[j]*y5[j]); z5[30*i+j] = sin(r)/r; } } /* plot 6 */ linspace(x6, -3, 3); linspace(y6, -4, 4); for(i=0; i<40; i++) { for(j=0; j<60; j++) { z6[60*i+j] = 3*(1-x6[i])*(1-x6[i])* exp(-(x6[i]*x6[i])-(y6[j]+1)*(y6[j]+1)) - 10*(x6[i]/5 - x6[i]*x6[i]*x6[i]-pow(y6[j],5))* exp(-x6[i]*x6[i]-y6[j]*y6[j]) - 1/3*exp(-(x6[i]+1)*(x6[i]+1)-y6[j]*y6[j]); } } subplot.subplot(2,3); /* create 2 x 3 subplot */ spl = subplot.getSubplot(0,0); /* get subplot (0,0) */ spl->title("Line"); spl->label(PLOT_AXIS_X,"t"); spl->label(PLOT_AXIS_Y,"Bessel functions"); spl->data2D(t, y1); spl->legend("j0", 0); spl->legend("j1", 1); spl->legend("j2", 2); spl->legend("j3", 3); spl = subplot.getSubplot(0,1); /* get subplot (0,1) */ spl->title("Polar"); spl->axisRange(PLOT_AXIS_XY, -1, 1, .5); spl->data2D(theta2, r2); spl->polarPlot(PLOT_ANGLE_RAD); spl->sizeRatio(-1); spl = subplot.getSubplot(0,2); /* get subplot (0,2) */ spl->title("3D curve"); spl->data3D(x3, y3, z3); spl->axisRange(PLOT_AXIS_X, 0, 400, 200); spl->axisRange(PLOT_AXIS_Y, -1, 1, 1); spl->axisRange(PLOT_AXIS_Z, -1, 1, 1); spl = subplot.getSubplot(1,0); /* get subplot (1,0) */ spl->title("Cylindrical"); spl->data3D(theta4, z4, r4); spl->coordSystem(PLOT_COORD_CYLINDRICAL, PLOT_ANGLE_DEG); spl->axisRange(PLOT_AXIS_Z, 0, 8, 2); spl->axisRange(PLOT_AXIS_XY, -4, 4, 2); spl->label(PLOT_AXIS_XYZ, NULL); spl = subplot.getSubplot(1,1); /* get subplot (1,1) */ spl->title("3D Mesh"); spl->axisRange(PLOT_AXIS_X, -10, 10, 5); spl->axisRange(PLOT_AXIS_Y, -10, 10, 5); spl->axisRange(PLOT_AXIS_Z, -.4, 1.2, .4); spl->data3D(x5, y5, z5); spl->colorBox(PLOT_OFF); spl->label(PLOT_AXIS_XYZ, NULL); spl = subplot.getSubplot(1,2); /* get subplot (1,2) */ spl->title("3D Mesh"); spl->data3D(x6, y6, z6); spl->axisRange(PLOT_AXIS_X, -3, 3, 2); spl->axisRange(PLOT_AXIS_Y, -4, 4, 2); spl->axisRange(PLOT_AXIS_Z, -8, 8, 4); spl->colorBox(PLOT_OFF); spl->label(PLOT_AXIS_XYZ, NULL); spl->ticsLevel(0.1); subplot.plotting(); }