// polygon.c --- n角形の1辺をN等分した場合の三角形分割のメッシュデータ生成
//   by 遠藤小欽
// 2015/04/20 作成
// 2015/05/15 改良版
// コンパイル
//   gcc -o polygon polygon.c
// 実行
//   ./polygon n N ファイル名

#include <stdio.h>
#include <stdlib.h>
#include <math.h>


#define Maxn (100) //正n角形
#define MaxN (250) // 一辺の長さをN等分
//f 正n角形の中心を０周，節点からなる正n角形の周数
//vx　正n角形の頂点のx座標
//vy　正n角形の頂点のy座標
//x　節点のx座標
//y　節点のy座標
//k　正n角形を8つの二等辺三角形に分け，その三角形の番号
//l　kブロックの何番目の節点番号または要素番号

// 実行　./polygon nの値　Nの値 メッシュファイルの名前

double x[Maxn*MaxN*(MaxN+1)/2+2],y[Maxn*MaxN*(MaxN+1)/2+2];
int b[MaxN+2],c[MaxN+2][Maxn];
double vx[Maxn+1],vy[Maxn+1];


int main(int argc, char **argv)
{
    int n,N,f,r,l,i,k,e;
    double theta,pi,t;
    FILE *fp;

    if (argc >= 3) {
      n=atoi(argv[1]);//atoi は文字列をintにする
      N=atoi(argv[2]);
    }
    else {
      printf("n: "); scanf("%d", &n);
      printf("N: "); scanf("%d", &N);
    }
    if (argc == 4) {
      fp = fopen(argv[3], "w");
      if (fp == NULL) {
        fprintf(stderr, "cannot open %s\n", argv[3]);
      }
    }
    else
      fp = stdout;

    fprintf(fp, "%d %d %d\n", 1+n*N*(N+1)/2, n*N*N, n*N);
    // 節点につけるラベルは、内点が0, 境界は1
    pi = 4.0 * atan(1.0);
    theta = pi/2;
    for (i=0; i<=n; i++) {
        vx[i]=cos(theta);
        vy[i]=sin(theta);
        theta -= 2*pi/n;
    }
    x[1] = 0;
    y[1] = 0;
    i = 2;
    fprintf(fp, "%f %f 0",x[1],y[1]);
    fprintf(fp, "\n");

    b[0] = 1;
    for(f = 1;f <= N;f++) {
        b[f] = 2 + (n*f*(f-1))/2;
        for(k=0;k<=n-1;k++) {
            for(l=0;l < f;l++) {
                e = 1 + k*(2*f-1) + l + n*(f-1)*(f-1);
                t = (double)l/f;
                x[i] = (1 - t)*(double)f/N*vx[k] +  t*(double)f/N*vx[k+1];
                y[i] = (1 - t)*(double)f/N*vy[k] +  t*(double)f/N*vy[k+1];
                if (f != N)
                    fprintf(fp, "%f %f 0\n",x[i],y[i]);
                else
                    fprintf(fp, "%f %f 1\n",x[i],y[i]);
                i++;
            }
        }
    }
    for (f=1;f<=N+1;f++) {
        b[f]=2 + (n*f*(f-1))/2;// f周目の最初の節点の節点番号
        for (k = 0; k<= n-1; k++) {
            c[f][k]=b[f]+k*f; // f周目の第kブロックの最初の節点番号
        }
    }

    for (k = 0;k < n-1;k++) {
        fprintf(fp, "1 %d %d 0\n",3+k,2+k);
    }
    fprintf(fp, "1 2 %d 0\n",n+1);
    /*
    fprintf(fp, "1 3 2 0\n"); fprintf(fp, "1 4 3 0\n"); fprintf(fp, "1 5 4 0\n");
    fprintf(fp, "1 6 5 0\n"); fprintf(fp, "1 7 6 0\n"); fprintf(fp, "1 8 7 0\n");
    fprintf(fp, "1 2 9 0\n");
    */
    for (f=2;f<=N;f++) {
        for (k=0; k<=n-1; k++) {
            for (l=1; l<=2*f-1;l++) {
                r=k*(2*f-1)+l;
                e=n*(f-1)*(f-1)+r;
                if (e == n*f*f) {
                    fprintf(fp, "%d %d %d 0\n", c[f-1][0],c[f][0],c[f+1][0]-1);
                }
                else if (e == n * f*f-1) {
                    fprintf(fp, "%d %d %d 0\n", c[f][0]-1,c[f-1][0],c[f+1][0]-1);
                    // c[f-1][0],c[f+1][0]-1,c[f][0]-1);
                }
                else if (l % 2 == 1) {
                    fprintf(fp, "%d %d %d 0\n",
                           c[f-1][k]+(l-1)/2,
                           c[f][k]+(l+1)/2, c[f][k]+(l-1)/2);
                }
                else { // l が偶数
                    fprintf(fp, "%d %d %d 0\n",
                           c[f-1][k]+(l-2)/2,
                           c[f-1][k]+l/2, c[f][k]+l/2);
                }
            }
        }
    }
    for (i=1; i<n*N; i++) {
        fprintf(fp, "%d %d 1\n", b[N]+i-1, b[N]+i);
    }
    fprintf(fp, "%d %d 1\n", b[N]+n*N-1, b[N]);
    fclose(fp);
    return 0;
}