10.3 別図

2011-09-30

10.3のFig.1と同様の図をプログラミング言語 Processingで描く。

Processingは、計算もできるし、図も描ける。8ビットパソコンに備わっていたBASICのように便利である。

プログラムの出力

出力結果をFig.2に示す。

プログラム

プログラムを以下に示す。

// _10_3graph3d.pde
// plot f(s,r) in 3D
// coded on 20110930

void setup() {
  size(600,450,P3D);
  translate(40,-40);
  xaxis();
  yaxis();
  zaxis();
  xnum();
  ynum();
  znum();
  noLoop();
}

float c = 20.0 / 12.0;
float fapr = 406.50, fjul = 416.64;
float fsep = 423.48, fdec = 433.84;

void draw(){
  float s, r;
  float [][] x   = new float[100][100];
  float [][] y   = new float[100][100];
  float [][] z   = new float[100][100];
  float [][] mx  = new float[100][100];
  float [][] my  = new float[100][100];
  float [][] mz  = new float[100][100];

  float sstart = 398.0;
  int       ns =  10;
  float rstart = 0.01;
  int       nr = 9;
  for(int i = 0; i < ns; i++){
    s = sstart + float(i);
    for(int j = 0; j < nr; j++){
      r = rstart + float(j)/100.0;
      x[i][j] = s;
      y[i][j] = f(s,r);
      z[i][j] = r;
      println("s r f=" + s + " " + f2(r) +" "+ y[i][j]);
    }
    println();
  }
  for(int i=0; i < ns; i++){
    for(int j=0; j < nr; j++){
      mx[i][j] = mapx(x[i][j]);
      my[i][j] = mapy(y[i][j]);
      mz[i][j] = mapz(z[i][j]);
      println("mx my mz=" + mx[i][j] +" "+ my[i][j] +" "+ mz[i][j] );
    }
  }
  for(int i=0; i < ns; i++){
    for(int j=0; j < nr-1; j++){
        line(mx[i][j],  my[i][j],  mz[i][j],
             mx[i][j+1],my[i][j+1],mz[i][j+1] );
    }
  }
  for(int j=0; j < nr; j++){
    for(int i=0; i < ns-1; i++){
      line(mx[i][j],  my[i][j],  mz[i][j],
           mx[i+1][j],my[i+1][j],mz[i+1][j] );
    }
  }
  save("10.3b.png");
}
// end of draw()

float mapx(float x){ return map(x, 398,    410,  0,      width); }
float mapy(float y){ return map(y,  -0,    700,  height, 0); }
float mapz(float z){ return map(z,   0.00, 0.08, 0,      -800); }

void xaxis(){
  line(mapx(398),mapy(0),mapz(0),
       mapx(408),mapy(0),mapz(0));
}
void yaxis(){
  line(mapx(398),mapy(0),  mapz(0),
       mapx(398),mapy(600),mapz(0));
}
void zaxis(){
  line(mapx(398),mapy(0),mapz(0.00),
       mapx(398),mapy(0),mapz(0.08));
}
void xnum(){
  fill(0);
  textSize(12);
  textAlign(CENTER,TOP);
  for(int i=390; i<409; i++){ text(str(i),mapx(i),mapy(0),mapz(0)); }
  text("S",mapx(409),mapy(0));
}
void ynum(){
  fill(0);
  textSize(14);
  textAlign(RIGHT);
  for(int i=0; i<600; i+=100){ text(str(i), mapx(398),mapy(i),mapz(0)); }
  text("r",mapx(398-0.25),mapy(550));
}
void znum(){
  textSize(14);
  for(int i=0; i<9; i++){
    float aa=i/100.0;
    text(str(aa), mapx(398),mapy(0),mapz(aa));
  }
}
float f(float s, float r){
  float[] a = new float[5];
  a[1]=f4(s,r)-fapr;
  a[2]=f7(s,r)-fjul;
  a[3]=f9(s,r)-fsep;
  a[4]=f12(s,r)-fdec;
  float sigma=0.0;
  for(int i=1; i<5; i++){ sigma+=a[i]*a[i]; }
  return sigma;
}
float f4(float s, float r){ return(s*(1+r/12)+c*(1+r/12)); }
float f7(float s, float r){ return(s*(1+5*r/12)+c*(5+15*r/12)); }
float f9(float s, float r){ return(s*(1+7*r/12)+c*(7+28*r/12)); }
float f12(float s,float r){ return(s*(1+10*r/12)+c*(10+55*r/12)); }
float f2(float x){ return ( round(x*100.0) / 100.0); }
float f3(float x){ return ( round(x*1000) / 1000); }

// end of program

history

2011-09-30 create.