Основная программа вызывается из файла cmas2d.bat и принимает в качестве параметра имя проекта. Это имя храниться в переменной projname.

Основная программа вызывает функции input (), calculate () и output ().

Функция input считывает файл .dat, сгенерированный GiD. Файл .dat содержит информацию о сетке. Функция 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");

Первая часть функции ввода соединяет имя проекта с расширением .dat, получая, таким образом, название файла данных. Этот файл открывается для чтения.

Далее вызывается функция jumpline(FILE*). Эта функция читает 80 байтов из файла, который она получает как параметр. Функция используется для перехода на следующую строку текста при чтении файла .dat.

for (i=0; i<6; i++) jumpline (fp);

fscanf(fp, "%d %d", &Nelem, &Nnod);

Первые шесть строк файла .dat пропускаются, так как это строки содержащие информацию для пользователя (см. .bas файл). После чего считывается общее количество элементов и узлов проекта и сохраняется в переменных 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);

В приведенном выше фрагменте кода выделяется память для хранения массива координат узлов (указатели x и y), массива связности (указатель N) и материалов, соответствующих каждому элементу (указатель imat).

В случае ошибки (нехватка памяти), создается файл с расширением .err. В этот файл выводится информация об ошибке и программа прерывается.

После чего следующие шесть строк пропускаются.

/* 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);

Далее считываются Координаты узлов и сохраняются в массивах x и y.

/* 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);

}

}

Информация о связности сохраняются в массиве N. Массив N логически можно представить в виде таблицы Nelem x 3. В массиве хранятся номера узлов (по 3 узла), для каждого элемента.

Для всех элементов считывается идентификатор материала. Если идентификатор материала - 0 (это означает, что для элемента материала не задан), создается файл .err, содержащий информацию об ошибке, и выполнение программы прерывается.

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);

}

Далее считывается оставшаяся информация из .dat файла.

Общее количество материалов сохраняется в переменной Nmat.

Плотность каждого материала сохраняется в таблице rho.

Общее количество условий считывается и сохраняется в переменной Ncnd.

Идентификаторы условий для узлов сохраняются в таблице nodc. Значения условий - в массиве wval.

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;

Эта функция вычисляет центр масс.

Также в нем приводится описание локальных переменных, используемых в calculate().

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);

Идентификаторы узлов текущего элемента сохраняются в переменных n1, n2, n3.

Для каждого элемента вычисляется объем (в данном случае вместо объема используется площадь, так как мы имеем дело с плоскими геометрическим объектом). Вычисленный объем сохраняется в переменной v.

Вычисляется геометрический центр элемента (совпадающий с центром масс). Его координаты сохраняются в переменных x_Cgi и y_Cgi.

Вычисляются суммы, которые войдут в числители координат центра масс. В конце результат деления переменных the x_num и y_num на переменную den сохраняются в переменные x_CG и y_CG.

void output() {

char filename[1024];

FILE *fp, *fplog;

double v;

Функция output() создает два файла: .post.res и .log.

Результаты для визуализации в постпроцессоре GiD сохраняются в файле .post.res. Именно этот файл хранит данные, которые позволяют GiD представить расстояние от соответствующего центра масс до каждой точки.

Вычисленное значение центра масс сохраняется в файле .log. Точность этого значения обратно пропорциональна размеру элемента.

/* 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);

Имя файла получается добавлением расширения .log к имени проекта, вычисленные значения координат центра масс, хранимые в переменных x_CG и y_CG, сохраняются в этом файле.

Cоздается файл .post.res. В этом файле хранятся результаты вычислений.

Формат файла .post.res file объясняется в справке GiD help, смотрите раздел

Posprocess 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);

В данном примере в файл .res записывается только скалярный результат с одним временным шагом.

Вот полный исходный код программы:

#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);

}