AUTO から,shooting point を得る.

wrtbv8(iap_type *iap, rap_type *rap, doublereal *par, integer *icp, doublereal *rldot, integer *ndxloc, doublereal **ups, doublereal **udotps, doublereal *tm, doublereal *dtm)
{
...
  mtot = ntot % 10000;
  fprintf(fp8,"%5ld",ibr);
  fprintf(fp8,"%5ld",mtot);
  fprintf(fp8,"%5ld",itp);
  fprintf(fp8,"%5ld",lab);
  fprintf(fp8,"%5ld",nfpr);
  fprintf(fp8,"%5ld",isw);
  fprintf(fp8,"%5ld",ntpl);
  fprintf(fp8,"%5ld",nar);
  fprintf(fp8,"%7ld",nrowpr);
  fprintf(fp8,"%5ld",ntst);
  fprintf(fp8,"%5ld",ncol);
  fprintf(fp8,"%5d\n",NPARX);

  // この下を追記.  
 // LABEL XXX の時のメッシュ点を splXXX.dat を出力する. 
  //========================================================
  // uniker : output the data of the orbit
  // ntst hooting points are got from AUTO
  // v[0],..., v[ntst-1]
  // v[ntst] sim v[0]
  char char_fname_sp[20];
  sprintf(char_fname_sp, "spl%03d.dat", mtot);
  FILE *fsp=fopen(char_fname_sp,"w");
  if(fsp == NULL) {
    printf("# ERROR: Could not open spXXX.dat");
    abort();
  }

  double period = par[10]; 
  //tm[j] in [0.0, 1.0]
  double s_j; // in [0.0, T]
  double h_j; // = s[j+1] - s[j]
  for (j = 0; j < ntst; ++j) {
    s_j = period*tm[j];
    h_j = period*(tm[j+1] - tm[j]); 
    fprintf(fsp,"    %ld", j);
    fprintf(fsp,"%19.10E", tm[j]);
    fprintf(fsp,"%19.10E", s_j);
    fprintf(fsp,"%19.10E", h_j);
    for (i = 0; i < ndim; ++i) {
      fprintf(fsp,"%19.10E", ups[j][i]);
    }
    fprintf(fsp,"\n");  
  }  

  s_j = period*tm[ntst]; // period*1.0
  h_j = period - s_j;  // 0.0 
  fprintf(fsp,"    %ld", ntst);
  fprintf(fsp,"%19.10E", tm[ntst]);
  fprintf(fsp,"%19.10E", s_j);
  fprintf(fsp,"%19.10E", h_j);
  for (i = 0; i < ndim; ++i) {
    fprintf(fsp,"%19.10E", ups[ntst][i]);
  }
  fprintf(fsp,"\n");  
  //========================================================
 // ここまで.この後は,元々の処理が続く
 
  /* Write the entire solution on unit 8 : */
  for (j = 0; j < ntst; ++j) {
    rn = 1. / ncol;
    for (i = 0; i < ncol; ++i) {
      k1 = i * ndim;
      k2 = (i + 1) * ndim;
      t = tm[j] + i * rn * dtm[j];
      fprintf(fp8,"    %19.10E",t);
      ...

なお,ups の中身は,以下のようになっている.

  // example
  // dimension = 2
  // 0<= j < ntst
  // (ncol-1) interpolated points
  //ups[j][0] // j_0 
  //ups[j][1]

  //ups[j][2] // j_1  
  //ups[j][3]
  //...
  //ups[j][(ncol-1)*dim]  // j_{ncol-1} 
  //ups[j][(ncol-1)*dim + 1] 
  
  //ups[j+1][0]
  //ups[j+1][1]
  //...
  //ups[ntst][0]
  //ups[ntst][1]