The main program

The main program is called from the cmas2d.win.bat file and has as parameter the name of the project. This name is stored in the variable projname.
The main program calls the input (), calculate () and output () functions.
The input function reads the .dat file generated by GiD. The .dat file contains information about the mesh. The calculate function read and processes the data and generates the results. The output function creates the results file.

void input () {
char filename[1024], fileerr[1024], sau1[1024], sau2[1024];
FILE *fp, *ferr;
int aux,j, error=0;
void jumpline (FILE*);
strcpy(filename, projname);
strcat(filename,".dat");
fp=fopen(filename,"r");

The first part of the input function concatenate the project name with the .dat extension, thus obtaining the name of the file that is to be read. This file is opened in order to be read.
The jumpline(FILE*) function is declared. This function simply reads a line from the file that it receives as a parameter, It is used to jump lines of the text when reading the .dat file.

for (i=0; i<6; i++) jumpline (fp);
fscanf(fp, "%d %d", &Nelem, &Nnod);

The first six lines of the .dat file are jumped over since these are lines of information for the user (see .bas file). Then the total number of elements and nodes of the project are read and stored in the variables Nelem and Nnod respectively.

x=(double *) malloc((Nnod+1)*sizeof(double)); if (x==NULL) {error=1;}
y=(double *) malloc((Nnod+1)*sizeof(double)); if (y==NULL) {error=1;}
N=(int *) malloc((Nelem+1)*3*sizeof(int)); if (N==NULL) {error=1;}
imat=(int *) malloc((Nelem+1)*sizeof(int));   if (N==NULL) {error=1;}
if (error) {
  strcpy(fileerr, projname);
  strcat(fileerr,".err");
  ferr = fopen(fileerr,"w");
  fprintf(ferr, "*****  ERROR: Not enough memory. *****\n");
  fprintf(ferr, "(Try to calculate with less elements)\n");
  fclose(ferr);
  exit(1);
}
for (i=0;  i<6; i++) jumpline (fp);

Space is reserved for storing the coordinates of the nodes (pointers x, y), the connectivities (pointer N), and the materials corresponding to each element (pointer imat).
In case of error (insufficient memory), a file is created with the extension .err. This file contains information about the error and the program is aborted.
The next six lines are jumped over.

/* reading connectivities */
for (ielem=1; ielem<=Nelem; ielem++){
  fscanf (fp, "%d", &aux);
  for(j=0;j<3;j++) fscanf (fp, "%d", &N[(ielem-1)*3+j]);
  fscanf (fp, "%d", &imat[ielem]);
  if (imat[ielem]==0){
    strcpy(fileerr, projname);
    strcat(fileerr,".err");
    ferr = fopen(fileerr,"w");
    fprintf(ferr, "**ERROR: Elements with no material!!**\n");
    fclose(ferr);
    exit(1);
  }
}

The connectivities are read and the N variable is saved. This variable is a Nelem x 3- size table with two fields. The nodes (assumed triangles of 3 nodes) forming the element are saved in the first field. The element identifiers are saved in the second one.
All the elements are checked, ensuring that they have been assigned a material. If the identifier of the material is 0 (meaning that no material has been assigned to the element), an .err file is created containing information about the error and the program is aborted.

for (i=0;  i<5; i++) jumpline (fp);
fscanf(fp, "%s %s  %d",sau1, sau2, &Nmat );
for (i=0;  i<3; i++) jumpline (fp);
/* reading density of each material */
for (i=1; i<=Nmat; i++)
fscanf (fp, "%d %lf", &aux, &rho[i]);
/* reading conditions*/
for (i=0;  i<4; i++) jumpline (fp);
fscanf(fp, "%d", &Ncnd);
for (i=0;  i<6; i++) jumpline (fp);
for (icnd=1; icnd<=Ncnd; icnd++) {
  fscanf (fp, "%d %lf", &nodc[icnd], &wval[icnd]);
  jumpline (fp);
}
fclose (fp);

Reading the remaining information in the .dat file.
The total number of materials is read and stored in the Nmat variable.
The density of each material are read and stored in the rho table. The material identifier indexes the densities.
The total number of conditions is read and stored in the Ncnd variable.
The nodes associated with a condition are read and stored in the nodc table indexed by the condition identifier. The value of the condition is stored in wval, another table indexed by the condition identifier.

void calculate ()
{
double v,aux1,aux2,aux3;
int n1, n2, n3;
int mat;
double x_CGi, y_CGi;
double x_num=0, y_num=0, den=0;
for(ielem=1;ielem<=Nelem;ielem++) {
  n1= N[0+(ielem-1)*3];
  n2= N[1+(ielem-1)*3];
  n3= N[2+(ielem-1)*3];
  /* Calculating the volume (volume is the area for surfaces) */
  v=fabs(x[n1]*y[n2]+x[n2]*y[n3]+x[n3]*y[n1]-x[n1]*y[n3]-x[n2]*y[n1]-x[n3]*y[n2])/2;
  x_CGi= (x[n1]+x[n2]+x[n3])/3;
  y_CGi= (y[n1]+y[n2]+y[n3])/3;
  mat= imat[ielem];
  x_num+= rho[mat]*v*x_CGi;
  y_num+= rho[mat]*v*y_CGi;
  den+= rho[mat]*v;
}
/* puntual weights */
for(icnd=1;icnd<=Ncnd;icnd++) {
  inod= nodc[icnd];
  x_num+= wval[icnd]*x[inod];
  y_num+= wval[icnd]*y[inod];
  den+= wval[icnd];
}
x_CG= (x_num/den);
y_CG= (y_num/den);

This is the function that calculates the center of mass.

The identifiers of the nodes of the present element are saved in n1, n2, n3.
This loop makes a rundown of all the elements in the mesh. The volume is calculated for each element. (Here, the volume is the area, provided we are dealing with 3D surfaces). The volume calculations are stored in the v variable.
The geometric center of the element is calculated (coinciding with the center of gravity) and the coordinates are stored in the x_Cgi and y_Cgi variables.
The numerator sums are calculated. When the loop is finished, the following sums are stored in the x_num and y_num variables. Finally, the result of dividing the x_num and y_num variables by the den variable is stored in the x_CG and y_CG variables.

void output() {
  char filename[1024];
  FILE *fp, *fplog;
  double v;

The output() function creates two files: .post.res, and .log.
The results to be visualized in GiD Post-process are stored in the .post.res file. It is this file that stores the data which enables GiD to represent the distance of each point from the corresponding center of mass.
The numerical value of the center of mass is saved in the .log file. The accuracy of this value is directly proportional to the element size.

/* writing log information file */
strcpy(filename, projname);
strcat(filename,".log");
fplog=fopen(filename,"w");
fprintf(fplog, "CMAS2D routine to calculate the mass center\n");
fprintf(fplog, "project: %s\n", projname);
fprintf(fplog, "mass center: %lf %lf\n", x_CG, y_CG);
fclose(fplog);

Creating the .log file: the .log extension is added to the project name and a file is created that will contain the numerical value of the position of the center of mass, which in turn is stored in the x_CG and y y_CG variables of the program.

Creating the .post.res file. The output data (results) are stored in this file.
The format of the .post.res file is explained in the GiD help, see section
Postprocess data files ->Postprocess results format.

/* writing .post.res */
strcpy(filename,projname);
strcat(filename,".post.res");
fp=fopen(filename,"w");
fprintf(fp,"GiD Post Results File 1.0\n");
fprintf(fp,"Result MC-DISTANCE \"LOAD ANALYSIS\" 1 Scalar OnNodes\n");
fprintf(fp,"ComponentNames MC-DISTANCE\n");
fprintf(fp,"Values\n");
for(inod=1;inod<=Nnod;inod++) {
  /* distance or each node to the center of masses */
  v=sqrt((x_CG-x[inod])(x_CG-x[inod])+(y_CG-y[inod])(y_CG-y[inod]));
  fprintf(fp,"%d %lf\n",inod,v);
}
fprintf(fp,"End values\n");
fclose(fp);

In this example only a scalar result , with a single time step, is written in the .res file.

This is the full source code of this program:

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

#define MAXMAT 1000
#define MAXCND 1000
char projname[1024];
int i,ielem,inod,icnd;
double* x,* y;
int* N,* imat;
int  nodc[MAXCND];
double rho[MAXMAT],wval[MAXCND];
int Nelem,Nnod,Nmat,Ncnd;
double x_CG,y_CG;

void input(void);
void calculate(void);
void output(void);

void main(int argc,char* argv[]) {
  strcpy(projname,argv[1]);
  input();
  calculate();
  output();
}

void input() {
  char filename[1024],fileerr[1024],sau1[1024],sau2[1024];
  FILE* fp,* ferr;
  int aux,j,error=0;
  void jumpline(FILE*);
  strcpy(filename,projname);
  strcat(filename,".dat");
  fp=fopen(filename,"r");
  for(i=0; i<6; i++) jumpline(fp);
  fscanf(fp,"%d %d",&Nelem,&Nnod);
  x=(double*)malloc((Nnod+1)*sizeof(double)); if(x==NULL) { error=1; }
  y=(double*)malloc((Nnod+1)*sizeof(double)); if(y==NULL) { error=1; }
  N=(int*)malloc((Nelem+1)*3*sizeof(int)); if(N==NULL) { error=1; }
  imat=(int*)malloc((Nelem+1)*sizeof(int));   if(N==NULL) { error=1; }
  if(error) {
    strcpy(fileerr,projname);
    strcat(fileerr,".err");
    ferr=fopen(fileerr,"w");
    fprintf(ferr,"*****  ERROR: Not enough memory. *****\n");
    fprintf(ferr,"(Try to calculate with less elements)\n");
    fclose(ferr);
    exit(1);
  }
  for(i=0; i<6; i++) jumpline(fp);
  /* reading the coordinates */
  for(inod=1; inod<=Nnod; inod++)
    fscanf(fp,"%d %lf %lf",&aux,&x[inod],&y[inod]);
  for(i=0; i<6; i++) jumpline(fp);
  /* reading connectivities */
  for(ielem=1; ielem<=Nelem; ielem++){
    fscanf(fp,"%d",&aux);
    for(j=0;j<3;j++) fscanf(fp,"%d",&N[(ielem-1)*3+j]);
    fscanf(fp,"%d",&imat[ielem]);
    if(imat[ielem]==0){
      strcpy(fileerr,projname);
      strcat(fileerr,".err");
      ferr=fopen(fileerr,"w");
      fprintf(ferr,"**ERROR: Elements with no material!!**\n");
      fclose(ferr);
      exit(1);
    }
  }
  for(i=0; i<5; i++) jumpline(fp);
  fscanf(fp,"%s %s  %d",sau1,sau2,&Nmat);
  for(i=0; i<3; i++) jumpline(fp);
  /* reading density of each material */
  for(i=1; i<=Nmat; i++)
    fscanf(fp,"%d %lf",&aux,&rho[i]);
  /* reading conditions*/
  for(i=0; i<4; i++) jumpline(fp);
  fscanf(fp,"%d",&Ncnd);
  for(i=0; i<6; i++) jumpline(fp);
  for(icnd=1; icnd<=Ncnd; icnd++) {
    fscanf(fp,"%d %lf",&nodc[icnd],&wval[icnd]);
    jumpline(fp);
  }
  fclose(fp);
}

void calculate() {
  double v;
  int n1,n2,n3;
  int mat;
  double x_CGi,y_CGi;
  double x_num=0,y_num=0,den=0;
  for(ielem=1;ielem<=Nelem;ielem++) {
    n1=N[0+(ielem-1)*3];
    n2=N[1+(ielem-1)*3];
    n3=N[2+(ielem-1)*3];
    /* Calculating the volume (volume is the area for surfaces) */
    v=fabs(x[n1]*y[n2]+x[n2]*y[n3]+x[n3]*y[n1]-x[n1]*y[n3]-x[n2]*y[n1]-x[n3]*y[n2])/2;
    x_CGi=(x[n1]+x[n2]+x[n3])/3;
    y_CGi=(y[n1]+y[n2]+y[n3])/3;
    mat=imat[ielem];
    x_num+=rho[mat]*v*x_CGi;
    y_num+=rho[mat]*v*y_CGi;
    den+=rho[mat]*v;
  }
  /* puntual weights */
  for(icnd=1;icnd<=Ncnd;icnd++) {
    inod=nodc[icnd];
    x_num+=wval[icnd]*x[inod];
    y_num+=wval[icnd]*y[inod];
    den+=wval[icnd];
  }
  x_CG=(x_num/den);
  y_CG=(y_num/den);
}

void output() {
  char filename[1024];
  FILE* fp,* fplog;
  double v;
  /* writing log information file */
  strcpy(filename,projname);
  strcat(filename,".log");
  fplog=fopen(filename,"w");
  fprintf(fplog,"CMAS2D routine to calculate the mass center\n");
  fprintf(fplog,"project: %s\n",projname);
  fprintf(fplog,"mass center: %lf %lf\n",x_CG,y_CG);
  fclose(fplog);
  /* writing .post.res */
  strcpy(filename,projname);
  strcat(filename,".post.res");
  fp=fopen(filename,"w");
  fprintf(fp,"GiD Post Results File 1.0\n");
  fprintf(fp,"Result MC-DISTANCE \"LOAD ANALYSIS\" 1 Scalar OnNodes\n");
  fprintf(fp,"ComponentNames MC-DISTANCE\n");
  fprintf(fp,"Values\n");
  for(inod=1;inod<=Nnod;inod++) {
    /* distance or each node to the center of masses */
    v=sqrt((x_CG-x[inod])*(x_CG-x[inod])+(y_CG-y[inod])*(y_CG-y[inod]));
    fprintf(fp,"%d %lf\n",inod,v);
  }
  fprintf(fp,"End values\n");
  fclose(fp);
  free(x);
  free(y);
  free(N);
  free(imat);
}

void jumpline(FILE* filep) {
  char buffer[1024];
  fgets(buffer,1024,filep);
}